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13.  ABSTRACT  (Maximum  2C0wcrcs> 

This  report  presents  the  results  of  several  studies  related  to  research  on  yield 
determination  and  discrimination  between  earthquakes  and  explosions.  The  work  pro¬ 
ceeded  along  two  lines:  (1)  the  determination  of  attenuation  of  Lg  coda  and  interme¬ 
diate-period  surface  waves,  and  (2)  discrimination  between  small  earthquakes  and 
explosions  in  the  central  United  States. 

The  attenuation  studies  in  the  first  four  chapters  of  this  report  conclude  that  (1) 
that  the  high  frequency  dependence  of  Lg  coda  Q  in  the  Basin  and  Range  can  be  caused 
by  a  rapid  Increase  of  Q  at  mid-crustal  depths,  (2)  that  patterns  of  upper  crustal  Q  ' 
and  upper  mantle  Q  correlate  throughout  the  western  Untied  States,  (3)  that  Lg  coda  Q 
values  in  the  southern  USSR,  China  and  southeastern  Asia  vary  between  200  and  700, 
and  (4)  that  attenuation  across  the  Arabian  Penninsula  is  high  for  a  shield  region.  ! 

Discrimination  studies  Indicate  (1)  that  a  neural  network  approach  can  be  effective 
in  determining  source  depth,  and  (2)  decreases  of  amplitude  with  frequency  in  the 
central  United  States  are  more  rapid  for  earthquakes  than  for  explosions. 
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Abstract 


Observed  Q  values  for  1-Hz  Lg  waves  which  traverse  the  Basin  and 
Range  province  exhibit  a  dependence  on  frequency  which  varies  as 
By  contrast,  the  frequency  dependence  of  intrinsic  Q  for  shear  waves  in  the 
crust  of  the  same  region  is  much  lower  and  may  be  nonexistent.  These 
apparently  conflicting  observations  can  be  explained  by  a  frequency- 
independent  layered  model  of  shear-wave  Q  which  has  very  low  values  in 
the  upper  crust  and  rapidly  increasing  values  at  greater  depths.  The  shear- 
wave  Q  model  which  best  explains  reported  values  of  QLg  and  its  frequency 
dependence  consists  of  an  8-km  thick  low-Q  layer  (60  or  less)  overlying  Q 
values  as  high  as  1000  or  more  at  mid-crustal  depths.  The  low  Qp  values  in 
the  upper  crust  are  best  explained  by  the  presence  of  fluids  in  interconnected 
cracks  and  pore  space.  At  greater  depths,  increasing  lithostatic  pressure 
closes  the  aacks,  leading  to  the  higher  Qp  values.  This  model  also  predicts 
fundamental-mode  Rayleigh  wave  attenuation  coefficients  which  agree  with 
observed  values.  The  result  of  this  study  contrasts  with  that  of  earlier  work 
in  high-Q  regions  of  the  eastern  United  States  where  shear-wave  Q  in  the 
upper  crust  must  vary  with  frequency  to  explain  the  observed  attenuation  of 
both  fundamental-mode  Rayleigh  waves  and  Lg. 

These  results  indicate  that  whereas  the  frequency-dependence  of  shear- 
wave  Q  is  an  intrinsic  property  of  crustal  material,  the  frequency  dependence 
of  QLg  arises  from  both  the  intrinsic  frequency  dependence  of  crustal  material 
and  its  layered  structure.  A  notation  is  suggested  which  distinguishes 
between  the  frequency  dependence  of  QLg  and  shear-wave  Q. 
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Introduction 


Data  that  can  be  used  to  study  the  attenuation  of  the  phase  Lg  and  its 
coda,  especially  at  frequencies  near  1  Hz,  has  increased  tremendously  in 
recent  years.  Observations  indicate  that  Q  for  these  waves  varies  regionally^ 
being  relatively  high  in  old  stable  regions  and  relatively  low  in  younger,  tec¬ 
tonically  active  regions  (Nuttli,  1973;  Singh  and  Herrmann,  1983;  Cong  and 
Mitchell,  1988).  Q  values  for  both  Lg  (QLg)  and  its  coda  (QJ  at  frequencies 
near  1  Hz  exhibit  similar  values,  being  within  10-15%  of  one  another  for  all 
regions  in  the  world  where  both  have  been  determined  (Nuttli,  1988).  The 
regional  variations  of  QLg  and  Qc  occur  in  patterns  which  vary  in  the  same 
way  as  regional  variations  in  crustal  shear  wave  Q  (Qp)  in  the  upper  crust 
obtained  from  surface  wave  attenuation  studies  (Mitchell,  1975;  Cheng  and 
Mitchell,  1981;  Hwang  and  Mitchell,  1987).  The  frequency  dependences  of 
QLg  and  Qc  also  vary  regionally,  but  in  the  opposite  sense  as  the  Q  values; 
they  are  lower  in  stable  regions  and  higher  in  tectonically  active  regions 
(NuttU,  1988). 

Since  Lg  is  a  guided  wave,  controlled  predominantly  by  crustal  proper¬ 
ties,  it  can  be  expected  that  the  regional  variation  of  QLg  and  Qc  and  their  fre¬ 
quency  dependences  are  caused  by  lateral  changes  in  the  anelastic  and 
scattering  properties  of  crustal  rock.  Models  of  Qp  which  explain  the  attenua¬ 
tion  of  both  Lg  waves  at  1  Hz  and  fundamental-mode  surface  waves  at  inter¬ 
mediate  frequencies  have  been  developed  for  several  regions  of  the  world, 
including  the  eastern  LFnited  States  (Mitchell,  1980)  and  western  United  States 
(Mitchell,  1981).  An  important  result  of  those  studies  is  that  the  frequency 
dependence  of  Qp  must  be  relatively  high  in  the  high-Q  crust  of  the  central 
and  eastern  United  States,  but  is  low  or  non-existent  in  the  low-Q  crust  of  the 
western  United  States  for  frequencies  of  about  1  Hz  and  lower.  Later  work 
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(Cong  and  Mitchell,  1988)  confirmed  the  earlier  results  in  North  America  and 
showed  that  the  same  conclusion  also  holds  for  South  America  and  India. 

The  results  described  above,  that  the  frequency  dependence  of  QLg  usu¬ 
ally  decreases  with  increasing  values  of  QLg  whereas  the  frequency  depen¬ 
dence  of  Qp  increases  with  increasing  Qp,  seem  to  be  contradictory.  The 
main  purpose  of  this  paper  is  to  explain  that  apparent  inconsistency.  It  will 
be  shown  that  the  frequency  dependence  of  QLg  in  the  Basin-and  Range  can 
be  produced  by  a  layered  Qp  structure.  Thus,  frequency  dependence  of  Qp, 
such  as  that  which  characterizes  crustal  rocks  of  the  eastern  United  States,  is 
not  required  in  the  Basin  and  Range.  A  standard  notation  will  be  proposed 
which  distinguishes  between  the  frequency  dependence  of  QLg,  Qc,  and  Qp. 
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Observations  of  Q^g  and  Qc  in  the 
Basin-and-Range 

Xie  and  Mitchell  (1990)  have  recently  determined  Q  and  its  frequency 
dependence  for  both  Lg  and  Lg  coda  waves  in  the  Basin-and-Range  province. 
Using  the  relation  QigiO  -  Qo^/  where  Qq  is  QLg  at  1  Hz  and  "n  is  a 
frequency-dependence  exponent,  they  found  that  QLg  can  be  described  by 
QLg(f)  =  (267  ±  56)f(0-37±0.06)  Lg  coda  Q  can  be  described  by  Qc(f)  = 

(275  ±  at  frequencies  between  0.2  and  2.5  Hz.  Their  determina¬ 

tions  utilized  a  stacked  spectral  ratio  (SSR)  method  (Xie  and  Nuttli,  1988) 
which  provides  stable  determinations  of  Q  with  variances  which  are  much 
smaller  than  those  associated  with  previously  used  methods.  Singh  and 
Herrmann  (1983)  obtained  Qc(f)  =  (250  ±  using  a  predominant 

frequency  analysis  and  Chavez  and  Priestley  (1986)  obtained  Qc(f)  =  (214  ± 

15)f(0.54i0.09) 

as  average  values  over  a  somewhat  broader  region  centered  in 
the  Basin-and-Range  province. 

The  similarity  of  the  Q  values  for  Lg  waves  and  coda  waves  and  the 
similarity  of  their  frequency  dependences  found  by  Xie  and  Mitchell  (1990) 
suggest  that  QLg  and  Q^.  can  be  used  interchangeably  when  using  those 
waves  to  infer  anelastic  properties  of  the  crust  or  in  using  the  attenuation  of 
one  or  the  other  of  those  waves  to  measure  magnitudes  of  earthquakes  or 
nuclear  explosions  from  Lg  waves,  at  least  in  low-Q  regions  like  the  Basin- 
and-Range  province.  That  similarity  is  particularly  convenient  when  deter¬ 
mining  crustal  models  of  Qp  which  explain  Lg  attenuation.  It  means  that 
amplitudes  of  Lg  waves  from  synthetic  seismograms  at  various  distances  can 
be  used  to  determine  QLg  values  predicted  by  crustal  models  and  these  can 
be  compared  with  observed  Q^.  values  which  are  more  abundant  than 
observed  QLg  values. 
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Implications  of  the  equivalence  of  QLg  and  were  discussed  by  Xie  and 
Mitchell  (1990).  They  noted  that  if  the  energy  flux  model  of  Frankel  and 
Wennerberg  (1987)  is  valid  at  distances  of  several  hundred  to  2000  km  in  the 
Basin  and  Range  province  the  equivalence  would  suggest  that  intrinsic  Q  is 
the  donunant  factor  affecting  wave  attenuation  there  and  the  effect  of  scatter¬ 
ing  Q  is  much  smaller.  Very  low  values  of  Qp  in  the  upper  crust  of  the 
Basin-and-Range  have  been  reported  by  Cheng  and  Mitchell  (1981),  Patton 
and  Taylor  (1984),  and  Lin  (1989)  from  studies  of  surface  waves  at  intermedi¬ 
ate  frequencies.  Because  of  the  relatively  long  wavelengths  used  in  the  sur¬ 
face  wave  studies  (20  -  150  km)  it  is  likely  that  those  low  Qp  values 
correspond  to  intrinsic  values.  Because  those  Qp  values  are  so  low  (50  -  100), 
it  is  likely  that  attenuation  over  the  frequency  range  0.2  to  2.5  Hz  used  by  Xie 
and  Mitchell  (1990)  is  also  dominated  by  intrinsic  Q  effects,  as  suggested  by 
Frankel  and  Wennerberg  (1987)  for  the  Anza  region  of  southern  California. 
In  the  following  sections,  we  consequently  assume  that  attenuation  at  both 
low  frequencies  and  frequencies  near  1  Hz  is  dominated  by  intrinsic  Q  and 
can  be  modeled  using  that  assumption. 
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Crustal  Qp  Models  for  the  Basin-and-Range 

Crustal  Qp  models  for  the  Basin-and-Range  province  have  recently  been 
obtained  by  Patton  and  Taylor  (1984)  and  Lin  (1989).  Although  they  both 
exhibit  Qp  values  which  are  low  compared  to  models  for  eastern  North  Amer¬ 
ica,  they  differ  from  one  another  in  the  distribution  of  those  values.  Qp“^ 
values  for  the  model  of  Patton  and  Taylor  (1984)  increase  continuously  from 
the  near  surface  to  upper  mantle  depths  whereas  the  model  of  Lin  (1989) 
includes  a  high-Q  zone  (low  Qp"^)  at  depths  between  about  12  and  25  km 
(Figure  1).  Both  models  were  obtained  using  the  assumption  that  Qp  is 
independent  of  frequency.  Frequency-independent  models  can  explain  the 
attenuation  of  both  fundamental-mode  surface  waves  and  1  Hz  Lg  waves  in 
that  region  (Mitchell,  1981).  Figure  1  also  shows  a  model  of  the  eastern 
United  States  (Cong  and  Mitchell,  1988)  which  must  be  frequency  dependent 
in  order  to  explain  the  attenuation  of  both  fundamental-mode  surface  waves 
and  Lg. 

Differences  between  the  two  Basin-and-Range  models  probably  occur 
because  the  two  studies  used  different  methods  to  determine  the  attenuation 
of  surface  waves.  Patton  and  Taylor  (1984)  used  the  method  of  Tsai  and  Aki 
(1969)  which  obtains  spectral  amplitudes  at  many  stations  surrounding  a 
seismic  event  and  solves  simultaneously  for  fundamental-mode  surface  wave 
attenuation  coefficients  and  the  moment  of  the  events.  The  method  is  most 
effective  in  regions  where  Q  has  little  lateral  variation,  but  is  susceptible  to 
systematic  errors,  which  can  give  attenuation  values  which  are  either  too  high 
or  too  low,  depending  on  path  configuration,  in  regions  where  Q  varies 
laterally  (Yacoub  and  Mitchell,  1977).  Observations  by  Singh  and  Herrmann 
(1983)  in  the  Basin-and-Range  indicate  that  Lg  coda  Q  varies  between  about 
200  and  300  in  that  region.  Thus,  the  model  of  Patton  and  Taylor  (1984) 
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could  be  biased  due  to  systematic  errors  in  surface  wave  attenuation  determi¬ 
nations. 

The  model  of  Lin  (1989)  was  obtained  by  inverting  two-station  measure¬ 
ments  of  both  fundamental-mo^’.e  and  first  higher-mode  Rayleigh  waves. 
Those  measurements  should  give  average  values  even  along  paths  where  Q 
varies  laterally.  That  model  includes  a  high-Q  zone  at  mid-crustal  depths. 
Hough  and  Anderson  (1988),  in  a  study  of  local  earthquake  recordings 
reported  a  similar  zone  of  high  Q  values  at  somewhat  shallower  depths  (5-12 
km)  beneath  the  Anza  array  in  southern  California. 

The  model  of  Lin  (1989)  will  be  used  as  a  starting  point  for  the  computa¬ 
tions  of  this  study.  Three  modified  forms  of  his  model  are  used,  aU  of  which 
include  a  low-Qp  layer  with  a  value  of  60  in  the  upper  crust.  Stacked  spectral 
ratios  are  computed  for  three  models,  one  in  which  the  low-Q  layer  is  6  km, 
another  8  km,  and  another  10  km  thick.  At  greater  depths  the  Qp  values  of 
Lin  (1989)  are  used. 
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Predicted  Values  for  QLg 

If  attenuation  in  the  Basin-and-Range  province  is  dominated  by  intrinsic 
effects,  then  the  attenuation  of  Lg  can  be  predicted  for  crustal  models  which 
have  been  determined  for  that  region.  In  the  present  study,  synthetic 
seismograms  for  the  Lg  phase  were  computed  at  several  distances  for  velocity 
and  Q3  models  for  the  crust  and  upper  mantle  of  the  Basin-and-Range  and 
the  attenuation  of  Lg  was  determined  using  the  stacked  spectral  ratio  (SSR) 
method  of  Xie  and  Nuttli  (1988);  predicted  Lg  Q  values  were  thus  determined 
using  the  same  method  that  was  used  to  obtain  observed  Q  values.  The  SSR 
method  was  applied  to  synthetic  seismograms  of  Lg  computed  by  mode  sum¬ 
mation  (e.g.,  Wang,  1981).  20  seismograms  were  computed  over  the  distance 
range  300  -  1250  km  for  Qp  models  and  velocity  models  obtained  by  Patton 
and  Taylor  (1984)  and  Lin  (1989),  Example  seismograms  computed  at  dis¬ 
tances  between  400  km  and  1200  km  appear  in  Figure  2  for  the  velocity  and  Q 
models  of  Lin  (1989).  The  rapid  loss  of  high  frequencies  caused  by  low-Q 
material  in  the  upper  crust  is  readily  apparent. 

The  solid  lines  in  Figure  3  denote  stacked  spectral  ratios  predicted  by 
three  modified  forms  of  the  model  of  Lin  (1989)  in  which  a  low-Qp  surface 
layer  overlies  much  higher  Qp  values  at  greater  depths.  The  calculated 
values  for  those  models  are  compared  to  observed  values  of  Xie  and  Mitchell 
(1990)  over  the  frequency  range  0.3  -  3.0  Hz.  The  best-fitting  model  is  that 
with  an  8  km  thick  low-Q  zone;  the  SSR  for  this  model  is  described  by  the 
relation  QLg  =  (248  ±  The  model  with  a  6  km  thick  low-Q  layer 

produces  higher  Q^  and  ti  values  and  the  model  with  a  10  km  thick  low-Q 
layer  produces  lower  values.  Results  for  the  model  with  an  8  km  thick  layer 
compare  weU  with  the  relation  QLg  =  (267  ±  obtained  by  Xie  and 

Mitchell  (1990)  from  stacked  spectral  ratios  of  observed  Lg  waves.  The  model 
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also  predicts  attenuation  coefficient  values  of  the  fundamental  Rayleigh  mode 
at  intermediate  periods  (6  -  35  s)  which  agree  satisfactorily  with  the  observed 
value  of  Lin  (1989)  as  shown  in  Figure  4,  The  only  discrepancy  is  at  a  period 
of  6  s  where  the  predicted  value  is  slightly  outside  the  one  standard  deviation 
range  of  the  observations.  The  model  with  a  10  km  thick  low-Q  zone  predicts 
values  which  are  in  even  better  agreement  in  the  6  -  7  s  range,  whereas  the 
model  with  the  6  km  thick  layer  predicts  values  which  are  substantially  lower 
than  observed  values.  Because  it  provides  the  best  fit  to  both  Lg  and  funda¬ 
mental  Rayleigh  wave  data,  the  model  with  an  8  km  thick  low-Q  zone  is 
taken  as  the  best  model  for  this  region.  The  dramatic  increase  of  Qp  at  mid- 
crustal  depths  is  most  easily  explained  as  being  caused  by  the  closure  of 
cracks  due  to  increasing  lithostatic  pressure.  As  suggested  by  Mitchell  (1975, 
1980),  the  low  Qp  values  in  the  upper  crust  imply  the  presence  of  fluid-filled 
cracks  and  pore  space.  Reductions  in  Q  due  to  fluid  movement  through  the 
crust  have  been  described  quantitatively  by  Mavko  and  Nur  (1975)  and 
O'Connell  and  Budiansky  (1977).  Brace  and  Kohlstedt  (1980)  summarized  the 
results  of  Brace  (1972)  and  Brace  (1980)  indicating  that  relatively  high  crustal 
permeability  extends  to  8  km  and  could  extend  to  greater  depths.  The  depth 
to  the  high-Qp  zone  obtained  in  this  study  is  therefore  consistent  with  depths 
commonly  associated  with  enhanced  crustal  permeability.  The  dashed  line  in 
Figure  3  was  produced  by  the  model  of  Patton  and  Taylor  (1984).  The 
absence  of  a  high-Q  layer  in  that  model  causes  observed  QLg  values  and 
frequency-dependence  values  to  be  too  low. 

The  modified  Qp  model  of  Lin  (1989)  which  produced  excellent  agree¬ 
ment  with  observed  values  of  Qi^g  and  its  frequency  dependence  is  a 
fi:equency-independent  model.  This  result  indicates  that  the  firequency 
dependence  of  Qj^g  near  1  Hz  which  is  observed  in  the  Basin  and  Range  can 
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be  produced  by  a  layered  structure  in  which  Qp  increases  rapidly  with  depth 
at  mid-crustal  depths.  It  is  possible,  however,  that  Qp  in  the  deeper  portions 
of  the  model,  particularly  the  high-Q  region,  can  be  frequency  dependent. 
When  Qp  is  high,  resolution  of  Qp  values  is  relatively  poor,  so  that  frequency 
dependence  caimot  easily  be  detected  (Mitchell,  1980). 

A  high-Q  layer  was  also  inferred  by  Hough  and  Anderson  (1988)  using 
body  waves  from  local  earthquakes  at  frequencies  between  15  and  100  Hz  for 
the  Anza  region  in  southern  California.  Their  Qp  model  consisted  of  a  high- 
Q  layer  (>  2000)  at  depths  between  5  and  12  km  sandwiched  between  a  5-km 
thick  surface  layer  with  Qp  values  of  about  500  and  deeper  material  with  Qp 
values  of  about  800.  Those  values  are  much  higher  than  those  of  the  present 
study  even  though  Q,,  values  there  are  lower  than  those  of  the  Basin  and 
Range  (Singh  and  Herrmann,  1983).  Their  model  therefore  suggests  that  cru¬ 
stal  Qp  values  increase  with  frequency  beginning  at  frequencies  somewhere 
between  5  Hz  and  15  Hz. 

The  depth  range  of  the  high-Q  zone  in  the  model  of  Hough  and  Ander¬ 
son  (1988)  corresponds  to  the  range  over  which  most  earthquakes  occur  in 
that  region.  The  high-Q  zone  in  the  model  of  Lin  (1989)  lies  between  depths 
of  about  12  and  25  km.  The  depth  distribution  of  earthquakes  has  been 
determined  in  the  western  Great  Basin  in  a  region  near  the  western  extent  of 
our  region  of  study  (Vetter  and  Ryall,  1983).  The  greatest  number  of  earth¬ 
quakes  occur  at  depths  between  10  and  15  km  in  a  zone  which  overlaps  the 
high-Q  zone  of  the  model  of  Lin  (1989),  and  earthquakes  occur  in  smaller 
numbers  at  depths  down  to  the  base  of  the  high-Q  zone.  Thus,  it  is  possible 
that  the  depth  distribution  of  earthquakes  in  the  Basin-and-Range  province  is 
related  to  the  Qp  distribution  there. 
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A  Classification  for  the  Frequency  Dependence 
of  QLg.  Qc/  and  Qp 

Expressions  for  the  internal  friction  (Q~^)  of  surface  waves,  assuming  that 
damping  losses  are  small  enough  that  can  be  neglected,  were  given  by 
Anderson  et  al.  (1965).  The  equation  for  Rayleigh  waves  is 


■‘<“>',1  ((q 


Qit  *  7^ 


“1 


'  itl  UCr  ap,  (Cr  Set, 

where  the  subscript  1  is  the  layer  index,  the  subscripts  R,  a  and  p  identify 
Rayleigh,  compressional,  and  shear  waves,  respectively,  to  is  angular  fre¬ 
quency,  and  Cr  is  Rayleigh  wave  phase  velocity.  The  subscripts  to,  a,  and  P 
indicate  that  those  quantities  are  held  constant  when  the  partial  derivatives 
are  computed. 

This  equation  can  be  used  to  calculate  Qr"^  for  any  Rayleigh  mode  by 
using  the  appropriate  velocities,  partial  derivatives,  and  intrinsic  Q  values. 
Mitchell  (1980)  extended  this  expression  to  incorporate  shear-wave  Q  values 
which  vary  with  frequency  as 


QpKO  =  efi  (2) 

where  ^  is  an  exponent  which  describes  the  frequency  dependence  of  Qp. 
Nuttli  (1986)  assumed  that  the  frequency  dependence  described  by  Mitchell 
(1980)  for  Qp  also  pertained  to  QLg.  However,  equations  (1)  and  (2)  indicate 
that  frequency  dependence  of  Qr“^  can  be  caused  either  by  frequency- 
dependent  Qp  in  the  crust  or  by  an  appropriately  layered  structure  in  which 
Qp  is  independent  of  frequency.  Results  of  the  present  study  indicate  that 
the  layered  structure  in  the  Basin  and  Range  causes  Qlj  to  increase  with  fre¬ 
quency.  That  increase  occurs  because  higher  frequencies  include  a  greater 
number  of  higher  modes  than  low  frequencies  and  those  higher-order  modes 
sample  more  deeply  in  the  crust  than  lower  order  modes.  A  high-Q  lower 
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crust  will  therefore  produce  a  positive  frequency  dependence  of  QLg.  Thus, 
the  frequency  dependence  of  Qlj  and  Qp  should  be  described  by  two  dif¬ 
ferent  parameters.  I  therefore  propose  that  the  frequency  dependence  of  the 
attenuation  of  Lg  waves  be  described  by  the  parameter  tj  and  that  the  fre¬ 
quency  dependence  of  internal  friction  for  shear  waves  be  described  by  the 
parameter  C.  Qlj  would  then  be  described  by 

QLgCO  =  Qof"-  (3) 

If  the  crust  and  upper  mantle  were  comprised  of  flat  layers  with  laterally 
uniform  elastic  and  anelastic  properties,  equations  (1)  and  (2)  could  provide  a 
good  approximation  for  the  effect  of  anelasticity  on  the  attenuation  of  Lg 
waves  which,  as  indicated  by  equation  (3),  may  have  a  different  frequency 
dependence.  It  is  known,  however,  that  heterogeneities  distributed 
throughout  the  lithosphere  produce  a  long  coda  following  the  Lg  wave  which 
cannot  be  replicated  in  synthetic  seismograms  computed  for  plane-layered 
structures  (Aki  and  Chouet,  1975).  Since  the  Lg  coda  is  composed  of  scat¬ 
tered  waves  which  travel  longer  paths  than  that  taken  by  Lg,  it  appears  that 
the  frequency  dependence  of  the  attenuation  of  Lg  coda  should  be  described 
by  yet  a  third  parameter  which  I  tentatively  designate  ti'.  Then 

Qc(f)  =  Qo'f"’-  (4) 

In  low-Q  regions,  such  as  the  Basin-and-Range,  it  appears  that  very  low 

intrinsic  Q  material  in  the  crust  attenuates  scattered  waves  at  the  same  rate  as 
it  attenuates  the  direct  Lg  wave.  Thus  the  relations  Qo'  =  Qo  arid  ti'  =  t]  may 
be  appropriate  for  all  low-Q  regions.  Further  work,  using  methods  capable 
of  providing  stable  measurements  of  and  Qj.  and  their  frequency  depen¬ 
dences,  will  be  necessary  in  high-Q  regions  to  see  if  the  relations  Qq'  =  Qo 
and  Tj'  =  hold  everywhere. 


13 


Conclusions 


All  of  the  frequency  dependence  of  QLg  in  the  Basin  and  Range  province 
can  be  explained  by  a  layered  frequency-independent  model  in  which  Qp 
increases  rapidly  with  depth  from  very  low  values  (<  100)  in  the  upper  8  km 
of  the  crust  to  high  values  at  greater  depths.  It  is  possible  that  the  high-Q 
layer  is  frequency  dependent,  but  that  that  dependence  cannot  be  detected 
because  of  poor  resolution  in  high-Q  material.  The  change  from  low  to  high 
values  of  Qp  most  likely  occurs  because  of  the  closure  of  cracks  produced  by 
increasing  lithostatic  pressure  with  depth.  The  frequency  independence  of 
this  model  contrasts  with  results  from  the  central  and  eastern  United  States 
where  frequency  dependence  of  Qp  in  the  upper  crust  is  required  to  explain 
the  attenuation  of  both  fundamental-mode  Rayleigh  waves  and  Lg  waves. 

Frequency  dependence  of  QLg  and  its  coda  at  frequencies  near  1  Hz  can 
therefore  be  caused  by  at  least  two  factors,  by  layering  in  a  frequency- 
independent  Qp  model  or  by  an  intrinsic  frequency  dependence  of  Qp.  It  is 
likely  that  both  factors  will  prevail  in  many  regions.  The  frequency  depen¬ 
dence  (ti)  of  QLg  and  the  frequency  dependence  (C)  of  Qp  are  thus  fundamen¬ 
tally  different  parameters  and  can  differ  greatly  from  one  another.  The  fre¬ 
quency  dependence  (ti')  of  Q^  may,  in  general,  be  different  from  both  ii  and 
5,  but  in  the  Basin  and  Range,  at  least  it  appears  that  ti  =  t|',  at  least  near  1 
Hz.  Since  that  equivalence  occurs  because  of  the  rapid  attenuation  of  scat¬ 
tered  waves,  it  is  likely  that  ti  =  in  other  low-Q  regions.  Further  work  is 
necessary  to  ascertain  whether  or  not  this  is  also  true  of  high-Q  regions. 
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Figure  Captions 


Figure  1.  Models  of  Qp"^  obtained  from  surface  wave  attenuation.  The  solid 
lines  indicate  a  frequency-dependent  model  for  the  eastern  United  States 
for  which  Qp  varies  as  f°‘®  (Cong  and  Mitchell,  1988).  The  dashed  Une  is  a 
frequency-independent  model  for  the  Basin-and-Range  province.  Both 
models  explain  the  attenuation  of  fundamental-mode  surface  waves  and 
1-Hz  Lg  waves  in  the  regions  they  represent. 

Figure  2.  Synthetic  Lg  seismograms  computed  by  model  summation  (Wang, 
1981)  for  distances  between  400  and  1200  km.  The  seismograms  were 
computed  for  along  an  azimuth  of  45°  from  a  vertical  strike-slip  fault  at  a 
depth  of  5  km  and  oriented  in  a  north-south  direction.  Seismograms  were 
computed  at  20  distances  and  stacked  to  obtain  stacked  spectral  ratios 
using  the  method  of  Xie  and  Nuttli  (1988).  The  time  scale  is  reduced  at 
each  distance  R  using  a  reduction  velocity  of  5.0  km/s. 

Figure  3.  Stacked  spectral  ratios  (circles)  obtained  by  Xie  and  Mitchell  (1990) 
for  the  Basin-and-Range  province  over  frequencies  between  0.3  and  3.0 
hz.  A  least  squares  fit  to  this  data  gives  QLg  =  (267  ±  56)f^°’**°'^\  The 
solid  lines  represent  theoretically  predicted  values  for  modified  versions  of 
the  model  of  Lin  (1989)  in  which  the  upper  crust  is  assigned  a  Qp  value  of 
60  and  allowed  to  vary  in  thickness  between  6  and  10  km.  The  dashed 
line  represents  theoretically  predicted  values  for  the  model  of  Patton  and 
Taylor  (1984). 

Figure  4.  Attenuation  coefficient  values  observed  by  Lin  (1989)  in  the  Basin 
and  Range  province.  Vertical  bars  denote  one  standard  deviation.  Solid 
lines  indicate  values  predicted  by  the  Qp  model  of  the  lower  crust  and 
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upper  mantle  in  Figure  1  combined  with  an  upper  crust  with  a  Qp  of  60. 
The  thickness  of  that  upper  crustal  layer  varies,  being  6,  8,  and  10  km 
thick  for  three  computations. 
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Chapter  2 


upper  Mantle  Anelasticity  and  Tectonic 
Evolution  of  the  Western  United  States 
from  Surface  Wave  Attenuation 


by 
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ABSTRACT 


Rayleigh  wave  phase  and  group  velocities  and  attenuation  coefficients  in 
the  period  range  18  to  120  sec  have  been  deterntined  for  three  regions  along 
nine  paths  across  the  western  United  States  using  a  two-station  technique. 
The  attenuation  coefficients  were  found  to  increase  from  east  to  west  between 
the  Rocky  Mountains  and  the  Pacific  coast. 

Rayleigh  wave  group  and  phase  velocities  were  inverted,  using  a  dif¬ 
ferential  procedure,  to  obtain  shear  wave  velocity  models  for  three  regions. 
These  velocity  models  were  then  used  in  an  inversion  process  where  Qp“^  as 
a  function  of  depth  was  obtained  from  observations  of  Rayleigh  wave 
attenuation.  The  inversion  results  show  that  Qp  values  in  the  upper  mantle 
of  the  western  United  States  are  lowest  for  the  coastal  regions  and  western¬ 
most  Basin  and  Range,  highest  for  the  Rocky  Mountains  and  western  Great 
Plains.  Intermediate  Qp  values  appear  to  occur  in  the  upper  mantle  beneath 
the  Basin  and  Range  and  the  Columbia  and  Colorado  Plateaus,  although 
uncertainties  in  the  data  prevent  a  clear  separation  between  that  region  and 
coastal  regions.  Low  values  of  Qp  occur  in  the  upper  crust,  higher  values  in 
the  lower  crust,  and  highest  values  in  the  uppermost  20  km  of  the  mantle  of 
all  these  regions.  These  Overlie  an  upper  mantle  low-Q  zone  in  which  Qp 
values  decrease  from  east  to  west,  lowest  values  (15-20)  being  similar  to  those 
observed  above  descending  slabs  in  the  western  Pacific. 

The  regional  pattern  of  upper  mantle  shear  wave  Q  values  is  similar  to 
that  observed  for  upper  crustal  shear  wave  Q  values  inferred  from  surface 
wave  studies  and  for  Lg  coda  Q  values.  The  patterns  of  upper  crustal  and 
upper  mantle  Q  variations  correlate  with  the  temporal  sequence  of  tectonic 
and  magmatic  activity  in  the  western  United  States,  lowest  Q  values  in  both 
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depth  ranges  occurring  where  tectonic  activity  has  been  most  recent.  Partial 
melting  and/or  enhanced  dislocation  motion  in  the  mantle  and  variations  in 
the  volume  of  interstitial  hydrothermal  fluids  in  cracks  in  the  upper  crust 
may  both  be  ultimately  due  to  processes  which  occurred  during  plate  con¬ 
sumption  and  interaction  near  the  western  margin  of  the  United  States. 
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INTRODUCTION 


Several  early  studies  revealed  that  the  upper  mantle  of  the  western 
United  States  is  characterized  by  abnormally  low  seismic  velocities  (Johnson, 
1967;  Archambeau  et  al.,  1969;  Julian,  1970).  Attenuation  studies  using  both 
body  waves  (Solomon  and  Toksoz,  1970;  Der  et  al.,  1980;  Der  et  al.,  1982)  and 
surface  waves  (Lee  and  Solomon,  1975;  Patton  and  Taylor,  1984;  Chen,  1985; 
Hwang  and  Mitchell,  1987)  have  shown  that  the  upper  mantle  there  is  also 
more  absorptive  of  seismic  energy  than  the  upper  mantle  of  more  stable 
regions.  These  studies  have  presented  strong  evidence  for  both  low-velocity 
and  low-Q  zones  in  the  upper  mantle;  however,  there  is  still  considerable 
uncertainty  about  their  regional  and  depth  distributions.  Moreover,  we  have 
little  understanding  of  the  physical  mechanisms  which  produce  lateral  hetero¬ 
geneity  of  the  zones.  It  is  presumably  associated  with  the  tectonic  develop¬ 
ment  of  the  region  and,  if  so,  can  be  correlated  with  surface  tectonics  or  with 
heterogeneity  in  the  crust  which  may  be  easier  to  observe  than  that  at  mantle 
depths.  The  main  thrust  of  the  present  paper  will  be  to  determine  the  extent 
of  heterogeneity  of  shear-wave  Q  (Qp)  across  the  western  United  States  using 
Rayleigh  waves  and  to  see  if  it  relates  to  variations  in  Q  obtained  at  crustal 
depths.  The  nature  of  that  relationship  may  provide  insight  into  the  nature 
of  evolution  of  the  upper  mande  and  crust  of  that  region. 

Early  studies  of  lateral  variations  in  mande  velocities  in  the  western 
United  States  utilized  body  waves  from  distant  sources.  Solomon  and  Toksoz 
(1970)  found  significant  lateral  variations  in  attenuation  of  long  period  P  and 
S  waves  beneath  the  United  States,  attenuation  values  being  high  between 
the  Rocky  Mountains  and  the  Sierra  Nevada  -  Cascade  ranges.  Lower 
attenuation  was  seen  at  stations  along  the  Pacific  coast.  Der  et  al.  (1975) 
investigated  the  regional  distribution  of  anelastic  attenuation  of  short  period 
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P  and  S  waves  and  found  that  the  observed  regional  distribution  pattern 
shows  high  attenuation  in  the  western  United  States,  including  California.  In 
more  recent  work.  Lay  and  Wallace  (1988)  observed  whole-mantle  attenuation 
using  ScS  waves.  They  found  Qscs  =  95  ±  4,  beneath  the  northern  Basin  and 
Range,  along  with  very  slow  shear  velocities.  They  found  much  lower 
attenuation  along  the  Pacific  Borderlands,  particularly  in  the  Pacific 
Northwest,  where  Qgcs  =  344  ±  88;  shear  velocities,  however,  were  found  to 
be  comparable  to  those  in  the  Basin  and  Range.  Their  results  suggested  that 
strong  variations  in  Qs^-s  within  tectonic  North  America  are  not  closely  corre¬ 
lated  with  shear  velocity  variations. 

Studies  of  shear  wave  Q  (Qp)  in  the  upper  crust  of  the  same  region  using 
surface  waves  have  shown  that  the  upper  crust  attenuates  seismic  energy 
much  more  severely  than  does  the  upper  crust  in  stable  regions  (Mitchell, 
1975;  Cheng  and  Mitchell,  1981;  Patton  and  Taylor,  1984).  The  high  anelastic 
attenuation  produced  by  those  low  Q  values  in  the  crust  also  produces  lower 
than  normal  values  of  Q  for  regional  phases  such  as  Pg  and  Lg  (Sutton  et  al., 
1967;  Singh  and  Herrmann,  1983). 

In  this  study  we  determine  Rayleigh  wave  attenuation  coefficients  over 
several  paths  in  the  western  United  States,  including  the  first  such  deterrai- 
nations  for  the  western  margin  region,  and  use  that  information  to  study  the 
lateral  variation  of  the  anelastic  properties  of  the  crust  and  upper  mantle. 
Surface  waves  should  be  more  sensitive  to  lateral  changes  in  rock  properties 
than  body  waves  since  body  waves  sample  the  mantle  in  a  vertical  sense  and 
will  be  more  sensitive  to  the  vertical  variations  than  to  the  lateral  variations  of 
rock  properties.  Another  advantage  of  using  surface  waves  in  such  studies  is 
that  different  wavelengths  sample  different  depths,  and  by  inverting  these 
measurements  we  can  constrain  the  depth  extent  of  the  attenuation  hetero- 
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geneity  beneath  the  western  United  States.  This  cannot  be  done  reliably 
using  body  waves  since  they  measure  attenuation  along  whole-mantle  paths 
(e.g.  Qscs)  or  average  differences  in  attenuation  along  mande-ray  paths  (e.g. 
P  and  S  waves).  In  this  study  we  emphasize  upper  mantie  variations  of  Qp, 
but  show  that  low  values  of  Qp  in  the  mande  and  upper  crust  appear  to  be 
closely  related  to  one  another  and  both  can  be  explained  by  a  simple  evolu¬ 
tionary  model  which  is  consistent  with  our  knowledge  of  the  tectonic 
development  of  the  western  United  States. 

DATA  ACQUISITION 

We  have  used  a  two-station  method  to  determine  fundamental-mode 
Rayleigh-wave  velocities  and  attenuation  coefficients  over  several  paths  in  the 
western  United  States  (Figure  1).  Nine  stations  belonging  to  the  World-Wide 
Standard  Seismograph  Network  (WWSSN)  and  one  station  of  the  Canadian 
Seismograph  Network  (CSN)  provided  data  for  this  study.  Since  we  res¬ 
tricted  our  study  to  Rayleigh  waves  only  long-period  vertical  component  data 
were  utilized. 

We  searched  for  earthquakes  which  were  recorded  by  two  stations  for 
which  the  great  circle  path  between  stations  differed  from  that  between  the 
epicenter  and  stations  by  4  degrees  or  less.  We  found  26  events  which  met 
that  criterion  and  were  well  recorded  at  two  stations  (Table  1). 

We  divided  the  western  United  States  into  three  regions  and  assumed 
that  elastic  and  anelastic  properties  of  the  crust  and  upper  mantle  are  rela¬ 
tively  laterally  uniform  within  each  of  the  regions.  The  paths  between  any 
two  stations  were  selected  so  that  they  are  totally  or  nearly  totally  located 
within  a  single  region.  In  order  to  minimize  the  effect  of  lateral  hetero- 


30 


geneities  we  tried,  as  far  as  possible,  to  select  paths  which  were  parallel  or 
subparallel  to  the  strike  of  the  tectonic  features. 

The  path  JCT-EDM  covers  the  first  region  which  we  designate  as  the 
Eastern  Rockies  (ER).  This  region  includes  a  small  part  of  the  Great  Plains  in 
addition  to  the  eastern  Rockies  area.  The  second  region,  designated  the 
Intermountain  region  (IM),  is  covered  by  6  paths,  LON-ALQ,  LON-DUG, 
COR-DUG,  DUG-TUC,  TUC-LUB,  and  TUC-JCT.  This  region  includes  the 
Basin  and  Range  province,  the  Colorado  and  Columbia  plateaus.  The  third 
region,  which  we  call  the  Western  Margin  (WM),  is  covered  by  two  paths, 
COR-GSC  and  COR-BKS.  This  region  includes  the  westenunost  Basin  and 
Range,  the  Sierra  Nevada,  the  Cascades,  and  the  Pacific  coast  ranges. 

DISPERSION  AND  ATTENUATION  MEASUREMENTS 

Analog  records  for  all  stations  were  digitized  and  processed  using  the 
multiple-filter  method  (Dziewonski  et  al,  1969),  after  proper  correction  for 
geometrical  spreading  and  instrument  response,  to  provide  group  velocities 
and  amplitude  spectra  for  surface  wave  modes  of  interest.  The  group  delay 
for  the  fundamental  mode  on  each  seismogram  was  used  to  calculate  the 
phase  which,  in  turn,  was  used  to  calculate  an  initial  estimate  of  wave 
number.  This  estimate  was  then  used  as  input  to  a  frequency-variable  filter 
(FVF)  (Russell  et  al.,  1988)  for  refining  the  spectral  amplitudes  and  group 
velocities. 

This  technique  combines  the  beneficial  aspects  of  a  time-variable  filter 
and  a  phase-matched  filter.  It  is  similar  to  the  time-variable  filter  (Landisman 
et  al.,  1969)  except  that  the  phase  of  the  mode  of  interest  is  first  removed  to 
minimize  spectral  distortions  due  to  phase  fluctuations. 
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After  frequency-variable  filtering  of  each  seismogram,  we  applied  a 
frequency-domain  Wiener  filter  (Hwang  and  Mitchell,  1986)  to  determine  the 
interstation  Green's  function  (G)  between  the  station  pairs  according  to  the 
relation 


G>f  (1) 

where  H  is  the  crosscorrelation  function  matrix  between  the  input  signal  (at 
station  1)  and  output  signal  (at  station  2),  and  R  is  the  autocorrelation  matrix 
of  the  input  signal.  Frequency-domain  Wiener  filtering  guarantees  the  remo¬ 
val  of  random  noise  and  provides  more  stable  and  accurate  estimates  of 
group  and  phase  velocities  and  attenuation  coefficients  than  other  methods, 
e.g.  the  spectral  ratio  technique. 

Knowledge  of  the  Green's  function  allows  us  to  measure  interstation 
phase  and  group  velocities  and  attenuation  coefficients.  Interstation  phase 
velocities  were  calculated  from  the  phase  spectra  of  the  transfer  functions 
using  the  formula  (Taylor  and  Toksoz,  1982) 


'  ft,  +  (*(0  ±  N) '  ^  * 

where  Ax  is  the  interstation  distance  in  kilometers,  t^,  is  the  first  time  point  of 
the  Green's  function  and  <{)(f)  is  the  phase  of  the  Green's  function  in  cycles. 
Following  Hwang  and  Mitchell  (1986),  attenuation  coefficients  of  the  intersta¬ 
tion  medium  were  calculated  using 


7(r)  =  -In  G(r) 
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where  A^  and  A2  are  the  epicentral  distances  in  degrees  for  stations  1  and  2, 
respectively,  and  G(f)  is  the  interstation  Green's  function.  The  interstation 
group  velocities  were  calculated  using 


U(f) 
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(4) 


where  and  X2  are  epicentral  distances  in  kilometers  for  stations  1  and  2, 
respectively,  and  tj  and  t2  are  arrival  times. 

Several  possible  sources  of  error  will  contribute  to  the  uncertainty  in  the 
measurements  of  phase  and  group  velocities  and  attenuation  coefficients  of 
surface  waves.  These  have  been  discussed  in  detail  by  several  authors  (e.g., 
Forsyth,  1975;  Kijko  and  Mitchell,  1983).  Systematic  errors  may  be  caused  by 
interference  by  other  phases,  incomplete  separation  of  modes,  and  lateral 
refraction.  Horizontal  refraction  in  an  inhomogeneous  medium  can  cause 
multiple  arrivals  of  energy  (Capon,  1970,  1971),  beating  (Pilant  and  Knopoff, 
1964),  or  unusual  amplitude  variations  (McGarr,  1969)  which  are  produced  by 
focusing  and  defocusing.  If  multiple  arrivals  are  out  of  phase  and  spatially 
separated,  then  phase-matched  filtering  (PMF)  will  take  care  of  these  and  iso¬ 
late  the  fundamental  mode  from  other  arrivals.  If,  however,  these  multiple 
arrivals  are  in  phase,  then  the  PMF  cannot  discriminate  between  them.  If 
that  occurs  energy  may  be  enhanced,  by  focusing,  at  one  station  relative  to 
the  other,  sometimes  resulting  in  negative  or  unrealistically  high  attenuation 
coefficient  values. 

Russell  et  al.  (1988)  showed  that  using  the  PMF  and  FVF  wiU  guarantee 
the  removal  of  the  phase  of  the  mode  of  interest  and  will  minimize  spectral 
distortions  due  to  phase  fluctuation.  Also,  the  amount  of  spectral  biasing 
due  to  windowing  (especially  when  there  is  significant  curvature  in  the 
amplitude  spectrum)  is  controlled  by  using  the  correct  window,  a  curvature 
correction,  and  changing  time  windows  appropriate  for  the  period  of  interest. 

Capon  (1970)  and  Forsyth  (1975)  demonstrated  that  most  late  arrivals  of 
surface-wave  energy  could  be  attributed  to  refraction  at  continental  margins. 
We  expect  the  same  situation  in  the  present  study  because  many  of  the 
events  used  in  this  study  are  located  in  the  Pacific  Ocean  near  the  Central 
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and  South  American  coasts.  Other  sources  of  systematic  errors  are  misloca- 
tion  of  the  source  and  uncertainties  in  origin  time.  For  two-station  measure¬ 
ments,  however,  such  errors  will  cancel.  Random  errors  can  be  caused  by 
errors  in  digitization  and  scattering. 

INTERSTATION  VELOCITIES 

The  Eastern  Rockies  Region  (ER) 

Five  two-station  events  along  the  path  JCT-EDM  were  used  to  determine 
interstation  phase  and  group  velocities  for  the  ER  region,  spanning  a  period 
range  from  13  to  92  sec.  Table  2  and  Figure  2a  show  the  means  and  standard 
deviations  of  these  five  events  which  will  be  used  later  to  invert  for  a  velocity 
li'.odel  of  the  region.  The  group  velocities  at  shorter  periods  exhibit  greater 
scatter  than  the  phase  velocities  and  the  group  velocities  at  longer  periods. 
This  may  occur  because  velocities  at  shorter  periods  are  more  sensitive  to  the 
lateral  variation  of  shallower  structure  than  longer  periods. 

The  Intermountain  Region  (IM) 

Thirteen  events  along  6  paths,  COR-DUG,  LON-DUG,  LON-ALQ, 
DUG-TUC,  TUC-LUB,  and  TUC-JCT,  were  used  to  determine  interstation 
phase  and  group  velocities  in  the  IM  region  for  the  period  range  16  to  80  sec 
(Table  3  and  Figure  2b).  Velocities  at  long  periods  show  scatter  which  is  just 
as  large  as  that  at  short  periods  for  this  region.  Scatter  at  longer  periods 
might  suggest  that  the  lateral  variations  in  structure  extend  to  very  great 
depths.  Large  standard  deviations  of  the  data  are  anticipated  in  this  region 
since  it  is  a  vast  area  which  covers  more  than  one  physiographic  province. 
However,  the  uncertainties  in  the  data  do  not  mask  the  main  features  of  the 
dispersion  curves;  they  are  sufficiently  dear  so  that  we  can  consider  the  velo- 
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city  structure  derived  from  these  dispersion  data  to  represent  an  average 
model  for  the  entire  region. 

The  Western  Margin  Region  (WM) 

For  this  region  we  selected  paths  which  are  parallel  or  subparallel  to  the 
strike  of  the  mountain  ranges  in  order  to  minimize  the  effect  of  lateral  varia¬ 
tions  in  the  structure.  Seven  two-station  events  were  used  along  two  paths, 
BKS-COR  and  COR-GSC,  to  determine  interstation  phase  and  group  veloci¬ 
ties  for  the  period  range  13  to  120  sec.  Table  4  and  Figure  2c  show  the  means 
and  standard  deviations  for  these  seven  events. 

The  largest  standard  deviations  in  the  dispersion  data  for  the  three  dif¬ 
ferent  regions  occur  for  the  IM  region,  whereas  they  are  lowest  for  the  ER 
region  and  intermediate  for  the  WM  region.  This  comparison  reflects  the 
complexity  of  the  IM  region,  in  which  multipathing  and  scattering  in  surface 
wave  propagation  could  be  greater  than  that  in  the  other  two  regions. 


INVERSION  FOR  VELOCITY  STRUCTURE 

The  phase  and  group  velocity  values  determined  earlier  have  been  used 
to  invert  for  shear  velocity  structure.  Bloch  et  ah  (1969)  showed  that  group 
velocities  are  more  sensitive  than  phase  velocities  to  small  differences  in  phy¬ 
sical  parameters.  Thus  inverting  group  velocities  can  produce  a  more 
detailed  velocity  structure  than  inverting  the  phase  velocities.  However,  it  is 
important  to  acquire  phase  velocity  information  because  group  velocities  are 
subject  to  an  additional  degree  of  nonuniqueness  compared  to  phase  veloci¬ 
ties  (Pilant  and  Knopoff,  1970).  We  have  therefore  simultaneously  inverted 
the  two  types  of  velocities  in  order  to  place  tighter  constraints  on  possible 
models  and  to  improve  structural  resolution  over  that  obtained  by  inverting 
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either  group  velocities  or  phase  velocities  alone  (Der  and  Landisman,  1972). 

In  the  present  study  we  have  adopted  a  half-space  as  the  initial  model  so 
that  our  final  model  is  not  biased  by  any  assumptions  made  concerning  the 
location  of  velocity  discontinuities  in  the  initial  model.  Numerical  experi¬ 
ments  show  that  stochastic  inversion  fails  to  give  a  stable  solution  when  start¬ 
ing  from  a  half-space  model,  but  the  differential  inversion  method  works 
well.  Thus  we  used  the  differential  inversion  approach  for  the  determination 
of  all  velocity  models.  Details  of  that  approach  are  given  by  Hwang  and 
Mitchell  (1987).  A  half-space  with  a  shear  velocity  of  5.0  knv/sec  was  used  as 
a  starting  model  for  aU  inversions.  Poisson's  ratio  was  taken  to  be  0.25  for 
the  crust  and  0.27  for  the  upper  mantle.  To  test  the  sensitivity  of  the  final 
velocity  model  to  initial  values  in  the  half-space,  we  also  ran  the  inversion 
process  with  starting  velocities  of  4.5  km/sec.  Velocities  in  the  final  model 
were  within  the  uncertainty  ranges  of  velocities  obtained  using  a  starting 
model  with  velocities  of  5.0  km/sec. 

Resolving  kernels  were  computed  for  all  models.  For  any  features 
derived  from  the  model  to  be  well  resolved,  the  width  of  the  resolving  ker¬ 
nels  should  be  narrower  than  the  thickness  of  that  feature  at  the  correspond¬ 
ing  depth. 

The  Eastern  Rockies  Region  (ER) 

Figure  3a  shows  the  velocity  model  along  with  its  resolving  kernels  for 
the  Eastern  Rockies  region  (ER).  The  theoretical  dispersion  values  predicted 
by  the  model  are  compared  to  observed  values  in  Figure  2a.  The  shear  velo¬ 
city  model  suggests  that  there  is  a  low-velocity  zone  in  the  upper  mantle 
which  seems  to  be  resolvable  since  the  width  of  the  resolving  kernel  is 
roughly  the  same  as  the  thickness  of  the  low-velocity  zone  at  the  correspond- 
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ing  depth.  This  low-velocity  zone  starts  at  a  depth  of  110  km  where  the 
shear  velocity  decreases  from  its  maximum  value  of  4.72  km/sec  and  reaches 
its  minimum  value  of  4.30  km/sec  at  250  km  depth.  There  is  about  a  9% 
decrease  in  shear  wave  velocity  compared  to  that  at  110  km  depth.  At  greater 
depths,  velocities  are  poorly  resolved,  since  resolution  decreases  with  increas¬ 
ing  depth. 

The  Intermountain  Region  (IM) 

Figure  3b  presents  the  shear-velocity  model  along  with  its  resolving  ker¬ 
nels  for  the  Intermountain  region  (IM).  The  theoretical  dispersion  values 
predicted  by  the  model  are  compared  to  observed  dispersion  values  in  Figure 
2b.  The  model  in  Figure  3b  has  large  standard  deviation  bars  throughout  all 
depths  because  the  large  standard  deviations  in  the  dispersion  data  have 
been  mapped  into  errors  in  the  model.  A  low-velocity  zone  starting  at  a 
depth  of  about  50  km  in  the  upper  mantle,  however,  seems  to  be  marginally 
resolvable  since  the  width  of  the  resolving  kernel  is  similar  to  the  thickness  of 
the  low-velocity  zone  at  the  corresponding  depth.  Shear  velocities  decrease 
from  4.5  km/sec  to  a  minimum  velocity  of  4.0  km/sec  at  a  depth  of  about  75 
km  and  continues  down  to  a  depth  of  about  130  km.  Because  of  the  poor 
resolution  at  those  depths,  it  is  possible  that  the  low-velocity  zone  could 
extend  considerably  deeper. 

The  Western  Margin  Region  (WM) 

For  the  Western  Margin,  the  upper  mantle  is  characterized  by  very  low 
velocities  of  about  4.1  km/sec  (Figure  3c).  Thus,  the  upper  mantle  lid  in  this 
region  may  be  either  very  thin  or  absent.  A  similar  result  was  found  by 
Archambeau  et  al.  (1969)  for  P  waves  in  the  Basin-and-Range  province  and 
Grand  and  Helmberger  (1984)  for  S  waves  through  a  broad  region  of  the 
western  United  States.  It  is  interesting  that  a  high-velocity  region  which 
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seems  marginally  resolvable  occurs  at  depths  between  180  and  250  km.  This 
might  be  a  high-velodty  lithospheric  plate  within  the  upper  mantle  which  is 
the  remnant  of  a  paleosubduction  zone.  The  presence  of  such  a  zone  in  the 
same  depth  range  has  been  identified  by  Koizume  et  al.  (1973). 

The  shear  velocity  models  in  Figure  3  indicate  no  systematic  regional 
changes  in  velocity  which  persist  at  all  depths,  but  at  depths  in  the  upper¬ 
most  mantle  (55-150  km)  where  resolution  is  relatively  good,  velocities  are 
higher  in  the  ER  region  than  they  are  in  the  two  other  regions.  The  similar¬ 
ity  in  average  velocities  for  the  IM  and  WM  regions  is  consistent  with  the 
results  of  Lay  and  Wallace  (1988)  who  found  that  both  Basin  and  Tvange  and 
Pacific  Borderland  have  comparable  velocities  and  are  associated  with  very 
slow  paths. 


RAYLEIGH  WAVE  ATTENUATION 

A  two-station  method  was  also  used  to  obtain  attenuation  coefficients  of 
fundamental-mode  Rayleigh  waves.  Interstation  attenuation  coefficients  were 
determined  from  the  interstation  Green's  function  as  described  earlier.  These 
measurements  were  much  more  difficult  than  those  for  phase  and  group 
velocities,  with  large  uncertainties  which  may  have  been  produced  by  focuss¬ 
ing  and  defocussing  effects  described  in  an  earlier  section.  Such  effects  have 
been  documented  by  Tanimoto  (1990)  for  part  of  the  region  of  this  study.  We 
attempted  to  minimize  this  problem  by  using  several  events  and  paths  for  the 
determination  of  attenuation  coefficients  and  by  processing  the  data  with 
modern  filtering  techniques.  Even  so,  some  of  the  data  were  characterized  by 
large  uncertainties  and  negative  attenuation  coefficient  values  sometimes 
occurred. 
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Focussing  of  waves  that  have  traveled  different  paths,  but  which  have 
nearly  the  same  path  from  epicenter  to  station,  results  in  an  increase  of 
amplitude  but  may  have  no  effect  on  the  wave  form.  In  addition,  because  of 
lateral  refraction,  the  spectral  amplitudes  observed  at  two  stations  may  not 
correspond  to  waves  which  left  the  same  portion  of  the  radiation  pattern  of 
the  earthquake.  If  this  occurs,  and  is  not  recognized,  amplitude  variations  at 
the  source  may  be  erroneously  attributed  to  th-i  effect  of  anelasticity  along  the 
path  between  stations.  For  these  reasons,  great  care  much  be  taken  in  select¬ 
ing  data  used  in  attenuation  studies.  One  method  of  checking  the  quality  of 
the  data  is  to  determine  the  particle  motion  of  the  waves  of  interest. 
Fundamental-mode  Rayleigh  waves  should  be  retrograde  elliptical  and  the 
plane  of  the  ellipse  should  be  along  the  propagation  path.  Plots  of  particle 
motion  were  made  in  several  cases  where  departures  from  a  great-circle  path 
were  suspected  and  data  were  discarded  if  observed  motion  departed  greatly 
from  that  expected  for  pure  Rayleigh  waves. 

The  Eastern  Rockies  Region  (ER) 

Attenuation  coefficient  values  of  the  fundamental  Rayleigh  mode  for  the 
ER  region  were  measured  for  five  events  at  periods  between  25  and  93  sec. 
The  means  and  standard  deviations  of  these  values  are  given  in  Table  5  and 
Figure  4.  Attenuation  values  at  short  periods  were  excluded  because  they 
exhibited  very  large  scatter.  For  periods  where  only  one  attenuation  coeffi¬ 
cient  value  is  available  we  assigned  the  maximum  standard  deviation  found 
for  other  periods  when  inverting  for  a  Q  model. 

The  Intermountain  Region  (IM) 

Interstation  attenuation  coefficient  values  for  the  IM  region  were  meas¬ 
ured  for  14  events  spanning  a  period  range  of  18  to  102  sec.  The  means  and 
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standard  deviations  of  these  values  are  given  in  Table  6  and  Figure  4.  The 
data  have  greater  scatter  in  the  period  range  18  to  50  sec  than  at  longer 
periods.  This  might  be  due  to  variations  in  crustal  structure  throughout  this 
region  since  paths  traverse  the  Basin  and  Range  province,  the  Columbia  Pla¬ 
teau,  and  the  Colorado  Plateau. 

The  Western  Margin  Region  (WM) 

Table  7  and  Figure  4  give  the  means  and  standard  deviations  of  the 
attenuation  coefficient  values  in  the  WM  region  measured  for  7  events  cover¬ 
ing  a  period  range  of  18  to  120  sec.  The  data  display  much  scatter,  which 
reflects  the  complexity  of  this  region,  especially  at  shorter  periods  which  sam¬ 
ple  the  crustal  structure. 

A  comparison  of  the  attenuation  values  in  the  three  regions  (Figure  4) 
shows  that,  without  considering  the  error  bars,  there  is  a  systematic  increase 
in  attenuation  coefficient  values  from  east  to  west  between  the  Rocky  Moun¬ 
tains  and  the  Pacific  coast.  If  we  take  the  uncertainties  of  the  estimates  into 
account,  even  though  there  is  much  scatter  in  the  data,  especially  at  short 
periods,  we  still  can  detect  some  systematic  variation  at  longer  periods.  This 
separation  is  obvious  between  the  ER  and  WM  regions  which  exhibit  the 
highest  and  lowest  values  throughout  the  long-period  range.  Values  for  the 
IM  region  lie  between  the  ER  and  WM  values  at  all  periods,  but  are  distinct 
from  them  at  only  a  few  periods. 

INVERSION  FOR  Qp"^ 

The  inverse  problem  for  obtaining  the  depth  dependence  of  Q  from  sur¬ 
face  waves  has  been  discussed  by  a  number  of  authors  such  as  Anderson  et 
al.  (1965),  Mitchell  (1975),  Lee  and  Solomon  (1975,  1978),  and  Hwang  and 


40 


Mitchell  (1987).  In  the  present  study  we  have  used  a  differential  inverse  pro¬ 
cedure  which  is  similar  to  the  one  used  for  obtaining  velocity  models  in  this 
study.  The  velocity  models  were  used  to  compute  partial  derivatives  of 
phase  velocity  with  respect  to  compressional  and  shear  velocities  which  are 
needed  for  the  inversion  process. 

An  initial  half-space  Qp  model  with  Qp  =  100  (Qp~^  =  0.01)  was  used  in 
the  inversion  process  for  all  regions.  This  value  is  based  on  results  of  a  pre¬ 
vious  study  in  the  Basin-and-Range  Province  (Patton  and  Taylor,  1984).  We 
tried  several  different  initial  Qp  values,  including  the  case  of  infinite  Qp  (Qp”^ 
=  0),  and  we  found  that  the  starting  model  had  only  a  small  effect  on  the 
resulting  Qp  models. 

The  inversion  process  was  applied  to  fundamental-mode  Rayleigh  wave 
attenuation  coefficients  values  chosen  from  Tables  5,  6  and  7  which  represent 
the  means  and  standard  deviations  of  the  ER,  IM,  and  WM  regions,  respec¬ 
tively.  To  avoid  anomalous  values  we  excluded  values  on  the  maxima  and 
minima  of  the  attenuation  coefficient  plots  and  selected  only  the  points  which 
we  considered  to  be  representative  of  the  attenuation  coefficient  mean  values 
for  each  region.  These  points  were  chosen  at  weU-separated  periods  so  that 
their  errors  could  be  considered  to  be  independent  (Der  et  al.,  1970).  There¬ 
fore,  for  the  inversion  of  attenuation  measurements  we  used  an  underdeter¬ 
mined  system  where  the  number  of  observations  is  less  than  the  number  of 
layers.  All  of  the  inversions  started  with  the  initial  Qp~^  value  of  0.01  and 
iterated  until  a  model  was  obtained  that  gave  the  best  fit  between  the  theoret¬ 
ical  attenuation  values  derived  from  that  model  and  the  observed  values. 

The  Eastern  Rockies  Region  (ER) 

The  Qp~^  model  obtained  for  the  ER  region  using  an  initial  Qp~^  value  of 
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0.01  is  shown  in  Figure  5  along  with  the  resolving  kernels  plotted  for  selected 
depths.  The  same  figure  also  shows  the  fit  between  the  theoretical  attenua¬ 
tion  coefficient  values  derived  from  the  resulting  model  and  the  observed 
values.  The  model  parameters  are  listed  in  Table  8.  For  upper  crustal 
depths,  Qp  values  are  about  70.  At  greater  crustal  depths  Qp  values  are 
higher,  ranging  between  100  and  250.  The  uncertainties  of  the  crustal  values 
are  likely  to  be  fairly  high,  both  because  of  the  lack  of  data  at  periods  below 
20  s  and  because  of  the  large  scatter  in  the  short-period  data.  Greater  confi¬ 
dence  can  be  assigned  to  the  upper  mantle  where  values  are  highest  (700  - 
850)  in  the  uppermost  20  km  at  depths  between  50  and  70  km.  A  low-Q  zone 
which  seems  resolvable  starts  at  a  depth  of  70  km  where  Qp  values  start  to 
decrease  and  attain  minimum  values  (~  35)  at  about  150  to  200  km  depths. 
At  depths  greater  than  200  km,  Qp  values  increase  slightly  but  resolution 
becomes  very  poor. 

The  model  resulting  from  the  inversion,  using  Qp~^  =  0  as  the  starting 
model  is  plotted  in  Figure  6,  along  with  its  resolving  kernels.  The  same  fig¬ 
ure  also  shows  the  fit  of  the  observed  attenuation  coefficient  values  with  the 
theoretical  values  predicted  from  the  model.  A  comparison  of  this  model 
with  that  in  Figure  5  shows  that  the  two  models  are  nearly  identical.  Thus 
the  differential  inversion  approach  assures  the  initial  Qp“^  values  will  not  bias 
the  resulting  Qp  model  even  if  those  values  are  zero. 

The  Intermountain  Region  (IM) 

Figure  7  shows  the  resulting  Qp“^  model  of  IM  region,  along  with  the 
resolving  kernels  plotted  for  selected  depths.  It  also  shows  the  fit  between 
observed  attenuation  coefficient  values  and  the  predicted  values  derived  from 
the  resulting  model.  The  large  standard  deviations  of  the  observed  data  are 
transformed  into  large  error  bars  in  the  resulting  model.  The  theoretical 
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attenuation  values  predicted  from  the  model,  however,  match  the  observed 
values  well.  The  model  parameters  are  listed  in  Table  9.  The  model  is 
characterized  by  low  Qp  values  (~  50)  at  upper  crustal  depths  and  relatively 
high  Qp  values  (~  130)  at  lower  crustal  depths.  As  for  the  ER  region,  the 
uncertainties  in  these  values  must  be  large.  The  highest  Qp  values  (~  700- 
300)  again  occur  in  the  uppermost  20  km  of  the  upper  mantle  at  depths 
between  30  and  50  km.  A  zone  of  very  low  Q  values,  which  seems  resolv¬ 
able,  starts  at  about  50  km  depth  and  attains  minimum  values  (===  20)  at  about 
175-250  km  depth.  At  greater  depths,  Qp  starts  to  increase  slightly  but  reso¬ 
lution  is  very  poor  at  those  depths. 

The  Western  Margin  Region  (WM) 

The  resulting  Qp~^  model  for  the  WM  region  is  characterized  by  very  low 
Qp  values  (~  30)  in  the  upper  crust,  and  higher,  but  still  low,  values  (~  90) 
in  the  lower  crust.  Again,  large  uncertainties  must  be  associated  with  those 
values.  The  highest  Qp  values  (about  300)  occur  in  the  uppermost  20  km  of 
the  upper  mantle  at  40-60  km  depth.  A  zone  of  extremely  low  Q  values, 
which  seems  resolvable,  starts  at  about  60  km  depth  and  attains  lowest  values 
(~  15)  at  about  175-250  km  depth.  At  larger  depths,  standard  deviations 
become  large  and  resolution  is  very  poor.  Table  10  lists  the  model  parame¬ 
ters  and  Figure  8  represents  the  model  with  its  resolving  kernels  plotted  for 
selected  depths.  Figure  8  also  shows  the  fit  between  the  observed  and  the 
theoretical  attenuation  coefficient  values  derived  from  this  model. 

A  comparison  of  the  Q3'’^  models  of  the  three  regions,  ER,  IM,  and  WM, 
is  presented  in  Figure  9  along  with  standard  deviations.  There  is  an  increase 
in  the  Qp“^  values  at  upper  mantle  depths  from  east  to  west;  Qp"^  values  are 
lowest  for  the  ER  region,  intermediate  for  the  IM  region  and  highest  for  the 
WM  region.  The  increase  in  Qp”^  values  between  the  IM  region  and  the  WM 
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region  at  upper  mantle  depths  is,  however,  within  the  overlap  of  the  error 
bars  and  cannot  strictly  be  concluded  to  occur.  The  difference  between  the 
IM  region  and  the  ER  region  is.  however,  unequivocal. 

DISCUSSION 

The  pattern  of  upper  mantle  values  of  Qp  obtained  from  Rayleigh  wave 
attenuation  coefficients  in  the  western  United  States  is,  in  a  gross  sense,  simi¬ 
lar  to  that  inferred  from  short-period  (1-Hz)  body-waves.  Der  et  al.  (1975) 
and  Der  et  al.  (1982)  found  upper  mantle  attenuation  in  the  southwestern 
United  States,  including  the  Basin-and-Range  and  California,  to  be  greater 
than  elsewhere  in  the  United  States.  The  small  number  of  available  stations 
and  the  wide  variation  in  individual  values,  however,  prevent  the  detection 
of  any  systematic  variations  between  the  Basin-and-Range  and  California. 
The  results  of  those  studies,  therefore,  while  consistent  with  our  results  that 
upper  mantle  Q  values  deaease  between  the  eastern  Rocky  Mountains  (ER) 
and  western  margin  (WM)  region,  neither  support  nor  argue  against  a 
decrease  in  upper  mantle  Q  values  between  the  Intermountain  region  (IM) 
and  western  margin  (WM). 

Solomon  and  Toksoz  (1970),  using  long-period  mantle  S  waves,  and  Lay 
and  Wallace  (1988)  using  long-period  ScS  waves,  found  upper  mantle  Q 
values  to  be  higher  to  the  west  of  the  Basin-and-Range  province  than  within 
it.  Lay  and  Wallace  (1988)  pointed  out  that  their  observed  pattern  of  attenua¬ 
tion  correlates  with  observed  heat  flow  in  the  western  United  States  and 
inferred  that  extensive  partial  melting  must  thus  occur  beneath  the  Basin- 
and-Range.  The  difference  in  the  inferred  pattern  of  upper  mantle  attenua¬ 
tion  obtained  using  long-period  body  waves  from  those  obtained  using 
short-period  body  waves  is  not  easily  explained,  but  is  probably  related  to  the 
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complex  structure  of  the  region.  One  or  both  results  could  be  biased  if  sta¬ 
tions  were  located  such  that  they  recorded  waves  which  traverse  local 
anomalous  features,  such  as  a  penetrating  slab  or  a  local  zone  of  partial  melt¬ 
ing.  The  number  of  stations  in  California  used  for  all  studies  was  small:  two 
for  direct  long-period  body  waves  (Solomon  and  Toksoz,  1970),  one  (and  one 
event)  for  long-period  ScS  (Lay  and  Wallace,  1988),  and  four  for  short-period 
body  waves  (Der  et  al.,  1982).  Resolution  of  this  question  may  be  possible 
when  sufficient  data  become  available  from  a  broad-band  network  currendy 
being  installed  in  southern  California. 

Our  results  from  inversions  of  Rayleigh  wave  attenuation  coefficient  data 
agree  better  with  patterns  of  attenuataion,  observed  for  short-period  body 
waves,  than  with  those  observed  for  long-period  body  waves.  The  results  of 
Der  et  al.  (1975)  and  Der  et  al.  (1982)  show  patterns  of  low  Q  associated  with 
the  upper  mande  over  a  broad  region  of  the  southwestern  United  States. 
Our  results  suggest  the  possibility  of  a  decrease  in  upper  mande  Q  values 
going  from  the  Basin-and-Range  province  westward,  but  there  is  much  over¬ 
lap  in  the  uncertainties  of  measured  surface  wave  attenuation  and  in  the 
resulting  models  of  Qp“^.  Consequendy  the  same  upper  mande  model  could 
satisfy  the  attenuation  data  (within  data  uncertainties)  in  both  the  IM  and 
WM  regions.  Changes  in  Qp~^  between  the  ER  and  IM  regions  are,  however, 
substantial  and  indicate  that  upper  mande  Q  values  must  decrease  more 
rapidly  between  the  eastern  Rockies  (ER  region)  and  the  Basin-and-Range  (IM 
region)  than  between  the  Basin-and-Range  and  the  western  margin  (WM) 
regions. 

Singh  and  Herrmann  (1983)  mapped  Q  values  for  1  Hz  Lg  in  the  United 
States  and  found  that  Q^g  decreases  from  about  800-850  to  about  200-300 
between  the  Rocky  Mountains  and  Basin-and-Range  and  decreases  further  to 
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about  150  in  coastal  California.  Thus,  crustal  values  of  shear  wave  Q  which 
control  Lg  Q  values  appear  to  decrease  drastically  between  the  Rocky  Moun¬ 
tains  and  Basin-and-Range,  but  decrease  much  less  severely  between  the 
Basin-and-Range  and  California.  This  pattern  resembles  the  variation  in 
upper  mantle  Q  values  inferred  from  our  study  of  surface  wave  attenuation. 
Evidence  for  decreasing  crustal  Q3  values  from  east  to  west  within  the 
western  United  States  were  first  found  by  Cheng  and  Mitchell  (1981)  in  a 
study  of  multi-mode  Rayleigh  waves.  Several  studies  using  various  types  of 
seismic  data  therefore  indicate  that  both  upper  mantle  Q  values  and  upper 
crustal  values  decrease  in  similar  patterns  from  east  to  west  in  the  western 
United  States.  The  similarity  in  these  pattern  suggests  that  the  anomalously 
low  values  of  Q  in  both  the  crust  and  upper  mantle,  as  well  as  the  systematic 
decrease  in  values  from  east  to  west,  have  a  single  underlying  cause. 

A  large  body  of  literature  has  linked  the  tectonic  regimes  and  magmatic 
activity  in  the  crust  and  mantle  of  western  United  States  to  plate  consump¬ 
tion  and  the  interaction  of  the  former  Farallon  plate  and  the  Pacific  plate  with 
the  North  American  plate.  The  Laramide  structure,  including  the  Rocky 
Mountains,  was  formed  generally  between  about  70  and  45  m.y.  B.P.  during 
the  early  Cretaceous  and  early  to  middle  Eocene  (Coney,  1972).  During  the 
mid-Cenozoic,  arc  magmatism  shifted  toward  the  Pacific  coast  as  the  angle  of 
plate  descent  beneath  the  Cordillera  changed  gradually  from  shallow  to  steep 
between  40  to  45  and  20  to  25  m.y.  B.P.  (Coney  and  Reynolds,  1977).  By  the 
beginning  of  Neogene  time,  a  continuous  magmatic  arc  extended  along  a  line 
parallel  to  the  continental  margin  in  the  western  Cordillera.  This  Neogene 
arc  was  later  disrupted  by  the  development  of  the  San  Andreas  transform 
(Dickinson  and  Snyder,  1978).  Along  the  Pacific  coast,  tectonic  activity  in  the 
Coast  Ranges  is  of  late  Cenozoic  age  (Hamilton  and  Myers,  1966)  and  both 


46 


tectonic  and  magmatic  activity  are  still  occurring  there. 

This  temporal  pattern  of  tectonic  and  magmatic  activity  corresponds  with 
the  spatial  pattern  of  Q  values  in  the  upper  crust  and  upper  mantle  of  the 
western  United  States.  The  correspondence  between  lowest  upper  mantle  Q 
values  and  recent  tectonic  and  magmatic  activity  suggests  that  those  regions 
are  sites  of  most  extensive  partial  melting  or  the  greatest  degree  of  dislocation 
motion  in  the  upper  mantle.  The  former  mechanism  is  favored  by  Shankland 
et  al.  (1981)  as  being  the  cause  for  low  seismic  velocities  and  Q  values  in  the 
upper  mantle,  whereas  the  latter  is  favored  by  Geuguen  and  Merrier  (1973) 
and  Minster  and  Anderson  (1980).  Because  different  attenuation  mechanisms 
are  additive  in  Q“^  (Shankland  et  al.,  1981),  it  is  possible  that  solid-state 
mechanisms  and  partial  nielt  mechanisms  act  together  to  produce  the  high 
attenuation  observed  in  the  upper  mantle  of  the  western  United  States. 

Mitchell  (1975,  1980)  suggested  that  Qp  in  the  upper  crust  is  controlled 
by  the  volume  of  interstitial  fluids  contained  in  cracks  and  microcracks.  Such 
fluids  could  be  hydrothermal  in  nature,  being  released  during  periods  of  tec¬ 
tonic  activity  when  temperatures  rose  to  high  levels  or  when  crustal  extension 
caused  reduced  pressure  at  crustal  depths.  If  that  is  the  case,  we  would 
expect  a  correlation  between  upper  crustal  Q  and  upper  mantle  Q  as 
observed  in  this  study.  It  is  interesting  that  those  two  low-Q  regions  are 
separated  by  a  zone  of  high  Q  values.  That  high-Q  region  may  occur  at 
depths  where  pressure  or  ductile  flow  closes  the  cracks  which  occur  at  shal¬ 
lower  depths  and  where  temperature  is  not  yet  high  enough  for  melting  or 
dislocation  motion  to  occur. 

At  upper  mantle  depths,  Qp  values  were  found  to  be  extremely  low, 
especially  for  the  IM  and  WM  regions  (minium  Qp  ~  15  -  20).  These  values 
are  comparable  with  other  upper  mantle  Qp  values  obtained  behind  several 
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island  arcs  in  the  Pacific  (Barazangi  and  Isacks,  1971;  Kanamori,  1970;  Moo¬ 
ney,  1970).  This  similarity  suggests  that  the  thermal  properties  and  tempera¬ 
ture  of  the  upper  mantle  in  western  North  America  remain  about  the  same  as 
they  were  when  active  subduction  was  occurring  during  the  Cenozoic  era. 


CONCLUSIONS 

Our  study  of  observed  surface  wave  dispersion  and  attenuation  in  the 

western  United  States  leads  to  the  following  conclusions: 

1.  Rayleigh  wave  attenuation  coefficients  and  upper  mantle  Qp“^  values 
increase  systematically  from  east  to  west  beneath  the  western  United 
States  between  the  Rocky  Mountains  and  the  Pacific  coast  in  a  pattern 
which  correlates  with  the  temporal  sequence  of  tectonic  episodes  and 
magmatic  activity.  Highest  upper  mantle  Qp"^  values  (lowest  Qp)  occur 
through  the  western  margin  of  the  United  States  where  the  occurrence  of 
tectonic  activity  is  most  recent,  whereas  lowest  values  (highest  Qp)  occur 
in  the  Rocky  Mountain  region  where  the  last  episode  of  tectonic  activity 
occurred  much  earlier. 

2.  The  attenuation  pattern  found  in  this  study  is  similar  to  that  of  short- 
period  teleseismic  body  wave  attenuation  patterns  determined  by  Der  et 
al.  (1975).  This  similarity  suggests  that  the  regional  variation  in  attenua¬ 
tion  of  both  short-period  body  waves  and  surface  waves  is  due  to  varia¬ 
tions  in  the  anelasticity  of  the  upper  200  km  of  the  mantle. 

3.  The  attenuation  pattern  is  also  similar  to  the  pattern  of  upper  crustal  Qp 
(Cheng  and  Mitchell,  1981)  and  of  Lg  coda  Q  at  about  1  Hz  (Singh  and 
Herrmann,  1983)  in  the  western  United  States.  This  suggests  that  the 
underlying  cause  for  Q  variations  in  the  upper  mantle  may  be  the  same 
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as  that  which  produces  lateral  variations  of  Q  in  the  upper  crust. 

4.  All  Qp  models  are  characterized  by  a  low-Q  zone  within  the  upper  man¬ 
tle  starting  at  about  50  -  70  km  depth.  Qp  values  at  those  depths  for  the 
Intermountain  and  Western  Margin  regions  are  extremely  low,  being 
similar  to  upper  mantle  Qp  values  obtained  behind  several  island  arcs  in 
the  western  Pacific. 
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TABLE  1.  Earthquakes  Used  in  This  Study 


Date 

Origin  Tune 
UT 

Latitude, 

deg 

Longitude, 

deg 

Depth, 

kin 

Ms 

Path 

April  20,  1964 

1156:41.6 

61.4 

-147.3 

30 

5.7 

COR-GSC 

July  11,  1964 

2025:40.3 

59.7 

-146.2 

40 

5.6 

COR-GSC 

March  28,  1965 

1322:57.6 

55.1 

162.1 

33 

5.9 

LON-DUG 

April  16,  1965 

2322:18.6 

64.7 

-160.1 

5 

5.8 

COR-GSC 

June  20,  1965 

1916:20.9 

25.4 

-109.4 

33 

5.8 

TUC-DUG 

Dec.  23,  1965 

2047:35.5 

60.6 

-140.7 

11 

5.8 

COR-GSC 

Aug.  18,  1966 

1033:17.7 

14.6 

-91.7 

85 

5.6 

ALQ-LON 

Nov.  3,  1966 

1624:32.7 

19.2 

-68.0 

33 

5.6 

JCT-TUC 

Nov.  4,  1967 

1626:48.2 

-2.8 

-77.7 

99 

6.0 

ALQ-LON 

Dec.  29,  1969 

051:47.2 

16.2 

-59.7 

17 

5.6 

JCT-TUC 

Nov.  18  1970 

2010:58.2 

-28.7 

-112.7 

33 

5.6 

BKS-COR 

April  18  1971 

1741:27.9 

-0.2 

-91.4 

33 

5.7 

JCT-EDM 

Aug.  20.  1971 

2136:9.6 

13.4 

-92.4 

33 

5.8 

JCr-EDM 

Aug.  2.  1972 

2138:50.4 

56.1 

163.2 

33 

5.6 

LON-DUG 

Aug.  12  1972 

1315:48.1 

5.0 

-82.6 

33 

5.7 

LON-DUG 

Oct.  20,  1972 

433:48.9 

20.6 

-29.7 

33 

5.7 

LUB-TUC 

June  22,  1974 

812:47.5 

-22.1 

-113.6 

33 

5.9 

BKS-COR 

June  27,  1974 

149:8.1 

33.8 

139.2 

16 

5.7 

DUG-COR 

Oct.  2  1974 

254:59.7 

-5.9 

-81.1 

5 

5.7 

LON-DUG 

Feb.  26,  1975 

448:54.7 

85.0 

98.0 

33 

5.4 

BKS-COR 

Aug.  19,  1975 

1457:12.3 

16.2 

-94.1 

85 

5.8 

JCT-EDM 

Nov.  21,  1975 

115:0.8 

7.5 

-77.6 

36 

5.8 

DUG-COR 

May  23,  1976 

1632:33.0 

-10.5 

-78.3 

73 

5.9 

LON-DUG 

April  4,  1978 

2111:41.8 

10.0 

-78.0 

35 

5.8 

DUG-COR 

June  4,  1979 

626:42.7 

15.7 

-93.6 

80 

5.7 

JCT-EDM 

July  21,  1979 

1212:14.0 

14.0 

-93.4 

26 

5.7 

JCT-EDM 
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TABLE  2.  Means  and  Standard  Deviations 
of  Fundamental  Mode  Rayleigh  Wave 
Velocities  in  the  Eastern 
Rockies  Region 


Period, 

s 

Group 

Velocity, 

km/s 

Standard 

Deviation, 

km/s 

Phase 

Velocity, 

km/s 

Standard 

Deviation, 

km/s 

13.0 

2.70 

0.56 

3.12 

0.08 

14.0 

2.73 

3.17 

0.04 

15.0 

2.71 

0.17 

3.20 

0.03 

16.0 

2.80 

0.31 

3.25 

0.01 

17.0 

2.71 

0.28 

3.27 

0.03 

18.0 

2.72 

0.27 

3.31 

0.05 

19.0 

2.71 

0.21 

3.36 

0.07 

20.0 

2.67 

0.16 

3.41 

0.08 

21.0 

2.58 

0.24 

3.46 

0.08 

22.0 

2.67 

0.19 

3.51 

0.08 

23.0 

2.88 

0.32 

3.56 

0.08 

24.0 

2.96 

0.29 

3.59 

0.08 

25.0 

2.98 

0.18 

3.63 

0.08 

26.0 

3.03 

0.17 

3.66 

0.08 

27.0 

3.09 

0.17 

3.68 

0.07 

28.0 

3.15 

0.17 

3.71 

0.07 

29.0 

3.19 

0.16 

3.73 

0.07 

30.0 

3.22 

0-15 

.  3.76 

0.08 

32.0 

3.27 

0.13 

3.80 

0.08 

34.0 

3.34 

0.07 

3.83 

0.08 

36.0 

4.42 

0.02 

3.86 

0.09 

38.0 

3.48 

0,04 

3.89 

0.09 

40.0 

3.54 

0.04 

3.91 

0.10 

43.0 

3.61 

0.04 

3.94 

0.10 

46.0 

3.69 

0.05 

3.96 

0.11 

48.0 

3.74 

0.05 

3.97 

0.12 

51.0 

3.83 

0.04 

3.98 

0.13 

53.0 

3.89 

0.04 

3.99 

0.13 

56.0 

3.98 

0.05 

3.99 

0.14 

59.0 

4.07 

0.08 

3.99 

0.16 

61.0 

4.05 

0.06 

4.07 

0.03 

64.0 

4.06 

0.05 

4.07 

0.04 

68.0 

4.06 

0.01 

4.07 

0.04 

72.0 

4.07 

0.01 

4.09 

0.04 

76.0 

4.08 

0.02 

4.08 

0.04 

80.0 

4.08 

0.02 

4.08 

0.05 

84.0 

4.07 

0.03 

4.08 

0.05 

88.0 

4.06 

0.04 

4.08 

0.06 

92.0 

4.05 

0.05 

4.08 

0.07 
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TABLE  3.  Means  and  Standard  Deviations  of 
Fundamental  Mode  Rayleigh  Wave  Velocities 
in  the  Intennountain  Region 


Group  Standard  Phase  Standard 

Period,  Velocity,  Deviation,  Velocity,  Deviation, 

s  km/s _ km/s _ km/s _ km/s 


16.0 

17.0 

18.0 

20.0 

21.0 

22.0 

23.0 

24.0 

25.0 

26.0 

27.0 

28.0 

29.0 

30.0 

32.0 

34.0 

36.0 

38.0 

40.0 

43.0 

46.0 

48.0 

51.0 

53.0 

56.0 

59.0 

61.0 

64.0 

68.0 

72.0 

76.0 

80.0 


3.05 

2.97 

3.03 

3.14 

3.09 

3.08 

3.08 

3.12 

3.17 

3.19 

3.22 

3.29 

3.37 

3.44 

3.54 

3.59 

3.61 
3.56 

3.60 

3.64 

3.62 

3.63 

3.65 
3.67 
3.69 

3.69 

3.70 
3.60 
3.73 
3.73 

3.76 

3.77 


0.04 

0.42 

0.38 

0.54 

0.41 

0.38 

0.35 

0.35 

0.35 

0.36 

0.37 

0.36 

0.37 

0.40 

0.45 

0.46 

0.41 

0.36 

0.38 

0.38 

0.37 

0.36 

0.36 

0.38 

0.41 

0.46 

0.52 

0.48 

0.40 

0.42 

0.50 

0.53 


3.30 

3.28 
3.32 
3.44 
3.47 
3.49 
3.52 
3.54 
3.57 
3.60 
3.62 
3.64 

3.66 

3.67 
3.70 
3.72 
3.74 
3.76 
3.78 
3.82 
3.87 
3.90 
3.94 
3.97 
4.01 
4.09 
4.14 

4.28 
4.07 
4.11 
4.20 
4.25 


0.10 

0.24 

0.21 

0.20 

0.20 

0.21 

0.21 

0.25 

0.24 

0.25 

0.25 

0.26 

0.27 

0.28 

0.28 

0.28 

0.29 

0.30 

0.31 

0.34 

0.39 

0.41 

0.44 

0.47 

0.52 

0.60 

0.69 

0.80 

0.54 

0.59 

0.72 

0.77 
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TABLE  4.  Means  and  Standard  Deviations 
of  Fundamental  Mode  Rayleigh  Wave 
Velocities  in  the  Western 
Margin  Region 


Period, 

s 

Group 
Velocity, 
km/s  - 

Standard 

Deviation, 

km/s 

Phase 

Velocity, 

km/s 

Standard 

Deviation, 

km/s 

13.0 

2.62 

0.04 

3.25 

0.03 

14.0 

2.87 

0.29 

3.31 

0.05 

15.0 

3.01 

0.43 

3.37 

0.10 

16.0 

2.91 

0.35 

3.38 

0.10 

17.0 

2.95 

0.22 

3.41 

0.12 

18.0 

2.95 

0.22 

3.45 

0.14 

19.0 

3.06 

0.25 

3.49 

0.14 

20.0 

3.12 

0.24 

3.51 

0.14 

21.0 

3.08 

0.21 

3.54 

0.15 

22.0 

3.18 

0.33 

3.56 

0.16 

23.0 

3.30 

0.40 

3.58 

0.17 

24.0 

3.33 

0.35 

3.59 

0.18 

25.0 

3.34 

0.30 

3.60 

0.20 

26.0 

3.34 

0.23 

3.61 

0.21 

27.0 

3.32 

0.18 

3.63 

0.22 

28.0 

3.29 

0.16 

3.64 

0.23 

29.0 

3.27 

0.19 

3.65 

0.24 

30.0 

3.27 

0.21 

3.67 

0.24 

32.0 

3.29 

0.23 

3.72 

0.22 

34.0 

3.31 

0.24 

3.75 

0.21 

36.0 

3.30 

0.22 

3.79 

0.21 

38.0 

3.39 

0.19 

3.82 

0.21 

40.0 

3.41 

0.17 

3.85 

0.22 

43.0 

3.48 

0.15 

3.80 

0.12 

46.0 

3.53 

0.13 

3.83 

0.12 

48.0 

3.55 

0.12 

3.80 

0.08 

51.0 

3.60 

0.13 

3.82 

0.08 

53.0 

3.58 

0.13 

3.83 

0.09 

56.0 

3.59 

0.13 

3.84 

0.09 

59.0 

3.60 

0.13 

3.86 

0.09 

61.0 

3.60 

0.14 

3.86 

0.09 

64.0 

3.61 

0.14 

3.88 

0.09 

68.0 

3.62 

0.15 

3.90 

0.09 

72.0 

3.65 

0.15 

3.91 

0.09 

76.0 

3.69 

0.15 

3.93 

0.09 

80.0 

3.73 

0.15 

3.94 

0.09 

84.0 

3.79 

0.15 

3.95 

0.09 

88.0 

3.84 

0.15 

3.96 

0.10 

92.0 

3.90 

0.15 

3.% 

0.10 

96.0 

3.97 

0.16 

3.96 

0.11 

100.0 

4.03 

0.16 

3.96 

0.12 

105.0 

4.09 

0.09 

4.05 

0.08 

110.0 

4.14 

0.09 

4.05 

0.09 
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TABLE  5.  Means  and  Standard 
Deviations  of  Attenuation 
Coefficients  in  the 
Eastern  Rockies 
Region 


Standard 

Deviation 


Period, 

s 

•yxlO-*, 

km"* 

xicr*, 

km"* 

25.3 

4.35 

1.75 

25.9 

4.40 

1.75 

26.6 

4.19 

1.75 

27.3 

3.85 

1.75 

28.1 

3.47 

1.75 

28.8 

3.10 

1.75 

29.3 

2.10 

0.72 

29.7 

2.77 

1.75 

30.1 

2.14 

0.42 

30.6 

2.48 

1.75 

31.0 

2.17 

0.29 

31.5 

2.26 

1.75 

32.0 

1.71 

1.05 

33.0 

1.81 

0.96 

34.1 

1.92 

0.85 

35.3 

1.62 

1.10 

36.6 

1.71 

1.03 

37.9 

1.76 

0.96 

39.4 

1.77 

0.88 

41.0 

1.72 

0.83 

42.7 

1.64 

0.88 

44.5 

2.01 

0.55 

45.5 

2.11 

1.75 

46.5 

2.11 

0.54 

47.6 

2.12 

1.75 

48.8 

2.19 

0.66 

50.0 

2.17 

1.75 

51.2 

2.17 

0.70 

53.9 

2.06 

0.63 

56.9 

1.71 

0.81 

60.2 

1.78 

0.86 

64.0 

2.05 

0.64 

66.1 

2.35 

1.75 

68.3 

1.81 

0.61 

70.6 

2.11 

1.75 

73.1 

1.87 

0.19 

75.9 

1.78 

1.75 

78.8 

1.65 

0.10 

81.9 

1.27 

1.75 

85.3 

1.30 

0.46 

89.0 

0.55 

1.75 

93.1 

0.22 

0.04 
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TABLE  6.  Means  and  Standard  Deviations  of  Attenuation 
Coefficients  in  the  Interxnountain  Region 


Period 

s 

■y  Standard  Deviation 

X 10-*  X 10-*  Period 

km  *'  km  s 

■y  Standard  Deviation 

X  10-*  X  10-4 

km  km  "* 

18.1 

10.49 

5.11 

33.0 

4.52 

2.79 

18.3 

10.13 

4.86 

33.6 

4.32 

2.99 

18.4 

9.68 

4.50 

34.1 

3.85 

2.79 

18.6 

9.11 

4.10 

34.7 

3.89 

2.77 

18.7 

8.42 

3.69 

35.3 

3.69 

2.33 

19.0 

7.61 

3.30 

35.9 

3.67 

2.14 

19.1 

6.68 

3.05 

36.6 

3.37 

1.78 

19.3 

5.63 

3.09 

37.2 

3.32 

1.81 

19.5 

5.58 

2.95 

38.6 

3.23 

2.34 

19.7 

5.52 

3.35 

39.4 

3.26 

2.45 

19.9 

5.07 

3.33 

40.2 

3.75 

2.86 

20.3 

4.21 

3.63 

41.0 

3.65 

3.01 

20.5 

4.59 

2.26 

41.8 

3.88 

3.45 

20.9 

4.96 

2.89 

42.7 

3.66 

3.50 

21.3 

4.86 

3.29 

43.6 

4.19 

3.90 

21.8 

4.44 

3.22 

44.5 

3.79 

3.64 

22.0 

4.83 

3.41 

45.5 

3.79 

3.70 

22.3 

3.84 

2.92 

46.5 

3.78 

3.28 

22.5 

3.94 

3.20 

47.6 

3.39 

3.08 

22.8 

3.16 

2.62 

48.8 

3.37 

2.67 

23.0 

3.00 

3.06 

49.9 

3.31 

2.52 

23.3 

2.82 

2.94 

51.2 

3.21 

1.96 

23.5 

2.93 

3.90 

52.5 

2.87 

1.72 

24.1 

2.58 

3.48 

53.9 

3.02 

1.43 

24.7 

2.31 

3.06 

55.4 

2.59 

1.13 

25.3 

2.15 

2.67 

56.9 

2.63 

1.22 

26.3 

2.95 

2.98 

58.5 

2.59 

1.47 

26.6 

3.38 

3.46 

60.2 

2.09 

1.09 

27.0 

3.36 

2.96 

62.1 

1.91 

1.22 

27.3 

3.83 

3.45 

64.0 

2.23 

1.19 

27.7 

3.77 

3.02 

66.1 

2.40 

1.17 

28.1 

4.28 

3.53 

68.3 

2.90 

1.18 

28.4 

4.17 

3.16 

70.6 

3.19 

1.36 

28.8 

4.74 

3.36 

73.1 

3.28 

1.76 

29.3 

4.56 

3.12 

75.9 

3.11 

2.19 

29.7 

4.52 

3.47 

78.8 

3.74 

2.58 

30.1 

4.47 

3.14 

81.9 

3.31 

2.23 

30.6 

4.82 

3.30 

85.3 

2.38 

1.41 

31.0 

4.69 

2.98 

89.0 

2.87 

0.01 

31.5 

4.68 

3.04 

93.1 

3.16 

1.31 

32.0 

4.53 

2.83 

97.5 

3.45 

1.81 

32.5 

4.73 

2.97 

102.4 

3.61 

1.55 
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TABLE  7.  Means  and  Standard  Deviations  of  Attenuation 
Coefficients  in  the  Western  Margin  Region 


Period 

s 

1 

X  10-* 

km  "* 

Standard  Deviation 

X  10-4 
km 

Period 

s 

y 

X  10-4 
km-> 

Standard  Deviation 
X  10-4 
km 

18.3 

12.21 

5.51 

30.6 

4.49 

3.54 

18.5 

10.60 

7.13 

34.1 

5.04 

3.98 

18.6 

11.00 

7.33 

34.7 

5.02 

3.99 

18.8 

11.36 

7.52 

35.3 

4.99 

3.99 

19.0 

11.60 

7.65 

35.9 

4.% 

3.99 

19.1 

11.67 

7.69 

36.6 

4.92 

3.98 

19.3 

11.56 

7.62 

37.2 

4.89 

3.96 

19.5 

11.34 

7.50 

37.9 

4.87 

3.92 

19.7 

9.50 

7.84 

38.6 

4.86 

3.85 

19.9 

9.28 

7.51 

39.4 

4.85 

3.74 

20.1 

9.02 

7.15 

40.2 

4.17 

3.49 

20.3 

8.73 

6.79 

41.0 

3.51 

3.21 

20.5 

8.42 

6.44 

41.8 

3.65 

3.00 

20.7 

8.09 

6.12 

42.7 

3.76 

2.82 

20.9 

7.74 

5.84 

43.6 

3.86 

2.68 

21.1 

7.40 

5.61 

44.5 

3.94 

2.60 

21.3 

7.05 

5.43 

45.5 

3.99 

2.57 

21.6 

6.71 

5.28 

46.5 

4.04 

2.61 

21.8 

6.38 

5.15 

47.6 

3.46 

2.71 

22.0 

6.06 

5.03 

48.8 

4.49 

2.12 

22.3 

5.75 

4.88 

50.0 

4.73 

1.93 

22.5 

6.40 

4.38 

51.2 

4.99 

1.76 

22.8 

6.10 

4.14 

52.5 

5.10 

1.98 

23.0 

5.80 

3.89 

53.9 

5.46 

1.97 

23.3 

5.49 

3.64 

55.4 

5.71 

1.99 

23.5 

5.18 

3.38 

56.9 

5.81 

2.03 

23.8 

4.86 

3.14 

58.5 

5.74 

2.08 

24.1 

4.54 

2.91 

60.2 

5.49 

2.09 

24.4 

4.24 

2.73 

62.1 

5.09 

1.99 

24.7 

3.95 

2.58 

64.0 

4.58 

1.73 

25.0 

4.28 

2.26 

66.1 

4.05 

1.29 

25.3 

4.08 

2.12 

68.3 

3.61 

0.69 

25.6 

3.89 

2.06 

70.6 

3.33 

0.11 

25.9 

3.70 

2.09 

73.1 

3.28 

0.55 

26.3 

4.17 

1.91 

75.9 

3.49 

0.94 

26.6 

4.02 

2.12 

78.8 

3.91 

1.22 

26.9 

3.86 

2.45 

81.9 

4.39 

1.53 

27.3 

3.50 

3.13 

85.3 

4.80 

1.97 

27.7 

3.62 

3.17 

89.0 

4.59 

2.30 

28.1 

3.74 

3.23 

93.1 

5.07 

2.67 

28.4 

3.86 

3.29 

97.5 

5.47 

3.15 

28.8 

3.98 

3.35 

102.4 

5.73 

3.57 

29.3 

4.10 

3.41 

107.8 

4.75 

6.06 

29.7 

4.22 

3.45 

113.8 

4.56 

5.67 

30.1 

4.35 

3.94 

120.5 

3.79 

4.23 
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TABLE  8.  Qj"*  Model  of  the 
Eastern  Rockies  Region 


Depth 

km 

Thickness 

km 

Q3 

Q3‘ 

X  10-^ 

Standard 

Deviation 

X  10^ 

5.0 

10.0 

65 

15.3 

3.0 

15.0 

10.0 

65 

15.3 

2.9 

25.0 

10.0 

74 

13.4 

2.6 

35.0 

10.0 

109 

9.1 

1.9 

45.0 

10.0 

226 

4.4 

1.3 

55.0 

10.0 

722 

1.4 

1.6 

65.0 

10.0 

834 

1.2 

2.3 

75.0 

10.0 

252 

4.0 

2.7 

85.0 

10.0 

111 

9.0 

2.8 

100.0 

20.0 

65 

15.3 

2.6 

120.0 

20.0 

46 

21.7 

2.1 

140.0 

20.0 

38 

26.4 

1.8 

162.5 

25.0 

35 

28.9 

2.1 

187.5 

25.0 

34 

29.1 

2.5 

212.5 

25.0 

36 

27.6 

2.6 

237.5 

25.0 

40 

24.9 

2.5 

262.5 

25.0 

46 

21.6 

2.2 

287.5 

25.0 

56 

17.9 

1.6 

312.5 

0.0 

71 

14.1 

1.0 

63 


TABLE  9.  Qp-' Model  of  the 
Intennountain  Region 


Depth 

km 

Thickness 

km 

Qp 

Qb' 

X  10-^ 

Standard 

Deviation 

X  10-^ 

5.0 

10.0 

45 

22.3 

8.2 

15.0 

10.0 

57 

17.7 

6.5 

25.0 

10.0 

133 

7.5 

3.7 

35.0 

10.0 

469 

2.1 

6.6 

45.0 

10.0 

334 

3.0 

9.8 

55.0 

10.0 

128 

7.8 

11.4 

65.0 

10.0 

71 

14.1 

11.4 

75.0 

10.0 

50 

20.1 

10.7 

85.0 

10.0 

39 

25.7 

10.0 

100.0 

20.0 

32 

31.0 

10.2 

120.0 

20.0 

27 

37.0 

11.0 

140.0 

20.0 

23 

43.6 

12.5 

162.5 

25.0 

20 

49.8 

14.2 

187.5 

25.0 

19 

54.1 

15.1 

212.5 

25.0 

18 

55.0 

14.9 

237.5 

25.0 

19 

52.3 

13.6 

262.5 

25.0 

22 

46.2 

11.4 

287.5 

25.0 

27 

37.4 

8.4 

312.5 

0.0 

38 

26.4 

4.9 

64 


TABLE  10.  Qp"*  Model  of  the  Western 
Margin  Region 


Depth 

km 

Thickness 

km 

Qs 

QS' 

X  10-^ 

Standard 

Deviation 

X  10-^ 

5.0 

10.0 

28 

35.7 

19.6 

15.0 

10.0 

29 

34.2 

18.2 

25.0 

10.0 

40 

25.1 

12.5 

35.0 

10.0 

86 

11.6 

9.0 

45.0 

10.0 

314 

3.2 

14.4 

55.0 

10.0 

297 

3.4 

18.2 

65.0 

10.0 

99 

10.1 

18.7 

75.0 

10.0 

50 

20.1 

17.4 

85.0 

10.0 

32 

31.2 

16.0 

100.0 

20.0 

24 

42.4 

15.5 

120.0 

20.0 

19 

52.4 

14.6 

140.0 

20.0 

16 

61.3 

15.4 

162.5 

25.0 

14 

68.8 

18.6 

187.5 

25.0 

14 

73.9 

21.8 

212.5 

25.0 

13 

75.8 

23.7 

237.5 

25.0 

14 

73.7 

23.7 

262.5 

25.0 

15 

67.6 

21.7 

287.5 

25.0 

17 

57.5 

17.9 

312.5 

0.0 

23 

43.7 

12.8 

65 


List  of  Figure  Captions 

Fig.  1.  Map  showing  two-station  paths  in  the  area  of 
study.  Triangles  denote  the  seismograph  stations. 

Fig.  2.  Interstation  fundamental  mode  Rayleigh 
phase  and  group  velocity  values  (squares)  for  the 
three  regions:  (a)  the  Eastern  Rockies,  (b)  the  Inter¬ 
mountain  region,  and  (c)  the  Western  Margin.  Verti¬ 
cal  bars  represent  standard  deviation.  The  solid 
curves  are  theoretical  dispersion  curves  predicted  by 
the  models  shown  in  Figure  3. 

Fig.  3.  Shear  velocity  model  (left)  obtained  by 
inverting  group  and  phase  velocities  of  the  funda¬ 
mental  Rayleigh  mode  for  the  three  regions:  (a)  the 
Eastern  Rockies,  (b)  the  Intermountain  region,  and 
(c)  the  Western  Margin.  Dashed  lines  indicate  the 
initial  half-space  model.  Resolving  kernels,  each 
normalized  to  unity,  corresponding  to  several  depths 
are  also  plotted. 

Fig.  4.  Comparison  of  attenuation  coefficients  of  the 
fundamental  Rayleigh  mode  for  the  three  regions: 
Western  Margin  region  (circles).  Intermountain 
region  (crosses),  and  Eastern  Rockies  region 
(squares).  Vertical  and  diagonal  lines  denote  stan¬ 
dard  deviation  values. 

Fig.  5.  model  (top  left)  obtained  by  inverting 

fundamental  Rayleigh  mode  attenuation  coefficients 
in  the  Eastern  Rockies  region.  Initial  =  0.01 

(dashed  line).  Normalized  resolving  kernels  at 
several  depths  are  also  shown.  Theoretical  attenua¬ 
tion  coefficients  (solid  line)  predicted  by  the 
model  are  compared  to  observed  values  (squares)  in 
the  lower  diagram.  Vertical  bars  represent  standard 
deviations. 

Fig.  6.  Qp~^  model  (top  left)  obtained  by  inverting 
fundamental  Rayleigh  mode  attenuation  coefficients 
in  the  Eastern  Rockies  region.  Initial  Qp"^  =  0.0. 
Normalized  resolving  kernels  (top  right)  at  several 
depths  are  also  shown.  Theoretical  attenuation  coef¬ 
ficients  (solid  line)  predicted  by  the  Qp~^  model  are 
compared  to  observed  values  (squares)  in  the  lower 
diagram.  Vertical  bars  represent  standard 
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deviations. 


Fig.  7.  model  (top  left)  obtained  by  inverting 

fundamental  Rayleigh  mode  attenuation  coefficients 
in  the  Intermountain  region.  Initial  Qp”^  =  0.01 
(dashed  line).  Normalized  resolving  kernels  (top 
right)  at  several  depths  are  also  shown.  Theoretical 
attenuation  coefficients  (solid  line)  predicted  by  the 
Qp~^  model  are  compared  to  observed  values 
(squares)  in  the  lower  diagram.  Vertical  bars 
represent  standard  deviations. 

Fig.  8.  Qp~^  model  (top  left)  obtained  by  inverting 
fundamental  Rayleigh  mode  attenuation  coefficients 
in  the  Western  Margin  region.  Initial  =  0.01 

(dashed  line).  Normalized  resolving  kernels  (top 
right)  at  several  depths  are  also  shown.  Theoretic^ 
attenuation  coefficients  (solid  line)  predicted  by  the 
Qp“^  model  are  compared  to  observed  values 
(squares)  in  the  lower  diagram.  Vertical  bars 
represent  standard  deviations. 

Fig.  9.  Comparison  of  Qp~^  models  of  Eastern  Rock¬ 
ies  (long-dashed  line).  Intermountain  (short-dashed 
line),  and  Western  Margin  (solid  line).  Horizontal 
bars  denote  one  standard  deviation. 
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Lg  Coda  Q  Variation 
Across  Eurasia 
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ABSTRACT 

Broad-band  digital  data  from  four  IRIS,  six  CDSN,  and  three  GDSN  sta¬ 
tions  in  the  USSR,  China,  and  southeast  Asia  have  been  collected  to  study  Lg 
coda  Q  in  Eurasia.  More  than  200  records  from  earthquakes  and  under¬ 
ground  nuclear  explosions  have  been  analyzed  with  a  stacked  spectral  ratio 
method.  Assuming  that  Lg  coda  Q  satisfies  the  relationship  Q(f)  =  Qo^/ 
where  Qq  and  “n  are  Q  at  1  Hz  and  the  power  law  fi-equency  dependence  of 
Q,  respectively.  A  tomographic  inversion  results  in  a  complex  pattern  of  Qq 
values.  Maximum  Qq  values  are  about  700,  relatively  low  for  a  stable  region. 
The  broadest  region  of  low  Qq  values  (200-300)  occurs  through  western 
China,  northern  Inndia  and  adjacent  regions  of  central  Asia  where  plate  colli¬ 
sion  has  produced  intense  tectonic  activity.  Other  regions  of  low  Qq  values 
lie  in,  and  east  of,  the  Caspian  Sea,  and  in  western  Turkey.  High  values  of 
Qo  (up  to  700)  occur  throughout  most  of  northern  Europe  and  Asia.  Bands 
of  moderately  high  values  of  Qq  (450-550)  lie  between  the  Black  and  Caspian 
Seas  and  in  eastern  China.  Frequency  dependence  values  show  no  clear  cur 
relationship  to  the  values.  They  are  highest  (1.0)  in  a  region  of  northern 
Siberia  where  Qo  takes  on  intermediate  values. 

Uncertainties  in  the  Qq  and  t]  estimates  are  relatively  high  in  some 
regions  such  as  the  Arabian  peninsula,  India,  eastern  China  and  the 
northwestern  portion  of  our  area  of  study.  It  is  expected  that  digital  stations 
currently  being  installed,  or  planned  for  the  near  future,  wiU  considerably 
improve  Qq  and  tj  estimates  for  some  of  those  regions  and  that  analog 
records  will  improve  them  for  others. 

INTRODUCTION 

Although  seismic  tomography  has  been  used  successfully  for  several 
years  to  study  lateral  variations  of  velocity,  it  has  only  recently  been  used  to 
study  lateral  variations  of  Q  (Xie  and  Mitchell,  1990).  This  development 
required  that  Lg  coda  Q  be  measured  with  relatively  high  precision.  Such 
precision  was  achieved  with  the  stacked  spectral  ratio  (SSR)  method  of  Xie 
and  Nuttli  (1988). 


Stacked  Spectral  Ratios 

The  SSR  method  is  applied  to  the  coda  of  Lg  waves  (Figure  1)  and  pro¬ 
vides  information  on  Q  averaged  over  the  areas  of  ellipses  (Figure  2)  which 
represent  areas  over  which  scattered  waves  sample.  We  used  a  total  of  more 
than  200  Lg  coda  time  series.  Each  of  the  coda  time  series  were  segmented 
by  many  windows  with  centered  lapse  times  starting  at  roughly  the  arrival  of 
group  velocity  of  3.15  km  s~^  and  increasing  to  successively  cover  the  entire 
signal  which  can  be  distinguished  from  the  ambient  noise.  We  used  a  10% 
Hanning  window  with  a  length  of  30  s  for  the  digitized  analogue  data  and  a 
length  of  25.6  s  for  the  digital  GDSN  data.  The  reference  noise  spectra  for 
each  record  were  obtained  in  the  way  described  by  Xie  and  Nuttli  (1988),  i.e., 
from  the  portion  of  the  trace  prior  to  the  P  arrival.  The  Fast  Fourier 
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transform  (FFT)  was  performed  on  both  signal  and  noise  and  the  effect  of 
noise  was  empirically  reduced  by  subtracting  the  reference  power  spectra  of 
noise  from  power  spectra  of  the  signal.  The  SSR  was  then  obtained  for 
single-trace  Q  measurements.  Figure  3  gives  some  examples  of  the  SSR 
analysis  over  the  frequency  band  from  0.2  to  2.0  hz. 

In  order  to  stablize  the  Q-inversion  using  the  random  coda  signal,  many 
windows  are  needed  in  the  SSR  method;  thus  we  must  use  the  later  part  of 
coda.  The  increased  stabUity  and  the  reduced  uncertainty  due  to  using  the 
later  coda,  however,  is  achieved  at  the  cost  of  poorer  spatial  resolving  power 
of  later  coda  since  it  samples  a  larger  area.  Eventually  one  achieves  an 
acceptable  uncertainly  in  coda  Q  inversion  with  limited  spatial  resolution. 

The  effect  of  random  noise  preserved  in  single-trace  Q  estimates  on  the 
^nal  image  can  be  empiorically  tested  utilizing  the  sample  standard  error  in 
Qn  due  to  the  randomness  of  SSRs  (Xie  ^d  Nuttli,  1988).  To  do  so  we 
denote  the  standaiJ  error  associated  with  by  SQ^,  n  =  1,  2,  3,  ...., 
(note  that  is  always  positive  due  to  the  way  standard  error  is  estimated), 
and  assume  that  SQ^  gives  a  good^ measure  of  the  absolute  value  of  real  error 
preserved  in  the  corrresponding  Q  measurements.  We  then  empirically  con¬ 
struct  a  number  of  noise  series,  whose  nth  member  has  an  absolute  value 
equal  to  and  a  sig^that  is  chosen  randomly.  The  nth  term  of  this  noise 
series  is  then  added  to  Q^,  which  we  shall  denote  as 

_  =  Qn  +  SQn 

Q'n  values  were  then  inverted  to  obtain  a  new  model.  The  difference 
between  the  new  image  of  Qp,  values  using  Q'^  and  the  original  image  of 
Qm  values  will  give  us  an  error  estimation  for  the  values.  This  process 
empirically  measures  the  effect  of  random  noise  on  the  image.  The  sign 
of  the  5Qn  series  can  be  simulated  by  pseudorandom  binary  generators  in  a 
computer  and  the  process  must  be  repeated  several  times  to  obtain  an  aver¬ 
age,  and  more  stable  error  measurement  of  values. 


Lg  Coda  Tomography 

Xie  and  Mitchell  (1990)  developed  a  back-projection  method  which  they 
applied  to  Lg  coda  to  obtain  a  tomographic  image  of  Lg  coda  Q  in  Africa. 
That  continent  was  divided  into  several  grids  (Ng)  where  was  assumed  to 
be  constant  within  the  nth  grid.  If  the  areal  average  of  Lg  coda  Q  (Q^)  in  the 
elliptical  area  sampled  by  coda  waves  is  obtained  at  maximum  lapse  time, 
then  denoting  the  area  of  the  nth  ellipse  that  overlaps  the  nth  grid  by  S^n  ^g, 
we  have 


1 

Qn 


+  E, 


n  = 


(1) 


where  is  the  summation  of  s^  over  the  ellipse  and  is  the  residual  due 
to  errors  in  modelling  Lg  coda  and  in  the  Q  measurement. 
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Lateral  Variations  of  Lg  Coda  Q 

The  back  projection  method  described  above  is  applied  to  image  lateral 
variations  of  Lg  coda  Q  in  central  and  southern  Eurasia.  We  have  divided 
the  whole  area  into  grids  with  widths  equal  to  in  latitude  and  longitude. 
The  grids  cover  about  75%  of  the  region  between  latitudes  10°  to  60°  north 
and  25°  to  135°  east. 

Figure  4  is  the  image  of  Lg  coda  Q  variation  at  1  Hz.  Several  features  are 
apparent  from  that  map.  A  broad  region  of  low  Qq  values,  between  200  and 
300,  spans  western  China,  northern  India,  and  adjacent  regions  of  central 
Asia.  More  localized  regions  of  low  Qq  values  lie  within,  and  near,  the 
Caspian  Sea  and  in  western  Turkey.  High  values  of  (up  to  700)  occur 
across  most  of  northern  Europe  and  Asia.  Bands  characterized  by  intermedi¬ 
ate  values  of  (450-550)  lie  between  the  Clack  and  Caspian  Seas  and  in 
eastern  China.  The  error  distribution  of  the  appears  in  Figure  5.  As 
expected  the  regions  with  largest  errors  occur  near  the  periphery  of  the  map 
where  path  coverage  is  sparsest.  Frequency  dependence  values  near  1  Hz 
were  determined  in  the  same  way  as  that  in  Africa  (Xie  and  Mitchell,  1990). 
Those  values  for  frequency  dependence  in  Figure  6  show  no  clear  relation¬ 
ship  to  the  Qo  values.  They  are  highest  (1.0)  in  a  portion  of  northern  Siberia 
where  Qq  values  are  high  and  lowest  (0.0)  in  southeast  Asia  where  Qq  takes 
on  intermediate  values.  It  should  be  noted,  however,  that  the  presently 
available  values  in  southeast  Asia,  as  well  as  some  others  (Figure  7),  are 
characterized  by  relatively  large  errors. 

Resolution  and  Errors 

The  imaging  of  Lg  coda  Q  variations  is  inherently  limited  because  of  the 
trade-off  between  the  stability  and  the  spatial  resolving  power.  For  this  rea¬ 
son  it  is  important  that  any  study  of  Lg  coda  Q  variations  include  an  analysis 
of  resolution  and  error.  We  have  calculated  point  spreading  functions  (psf) 
to  estimate  the  resolving  power  at  the  locations  where  these  grids  are  cen¬ 
tered.  All  of  the  psf  determinations  were  normalized  such  that  the  maximum 
values  are  unity.  In  general,  the  degree  of  spreading  of  psf's  increases  as  the 
density  of  coverage  decreases. 

Figure  8  shows  example  point  spreading  functions  centered  on  four 
widely  spaced  grids.  The  broadest  psf  occurs  in  the  northeastern  portion  of 
the  map  where  path  coverage  is  poor. 

The  psf  plots  allow  us  to  determine  whether  or  not  features  of  the  Qo 
and  iq  images  are  real  or  whether  they  are  simply  artifacts  of  poor  resolution. 
The  low-Q  region  just  south  of  the  Caspian  Sea  (Figure  4)  is  of  particular 
interest.  It  has  long  been  known  that  the  Caspian  Sea  is  characterized  by 
high  attenuation,  but  our  results  suggest  that  those  low  Q  values  extend  well 
to  the  south  of  the  sea.  The  psf  in  that  region  suggests  that  the  southern 
extension  of  that  low-Q  region  is  real  and  is  not  a  result  of  poor  resolution. 

Conclusions 
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Qo  values  throughout  central  and  southern  Asia  range  between  200  and 
700  and  frequency  dependence  values  range  between  0.3  and  nearly  1.0.  No 
clear  relation  is  apparent  between  those  values.  Lowest  values  occur  in 
and  south  of  the  Caspian  Sea,  near  the  Black  Sea,  and  in  a  broad  band  from 
India  to  western  China.  Highest  values  occur  in  northern  Asia. 

Uncertainties  in  Qq  and  q  estimates  are  high  in  regions  around  the  peri¬ 
phery  of  the  study  area.  Additional  data  with  better  station  coverage  will  be 
necessary  to  reduce  those  uncertainties. 
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Figure  Captions 


Figure  1.  Vertical  component  seismograms  from  three  stations  in  China.  The 
abscissa  indicates  the  lapse  time  measured  from  the  earthquake  origin 
time. 


Figure  2.  Sampling  pattern  of  Lg  coda  used  in  this  study.  Each  record  of  Lg 
coda  is  assumed  to  sample  an  elliptical  area  corresponding  to  the  max¬ 
imum  lapse  time  used  in  the  analysis.  When  the  continental  boundary  is 
encountered  it  is  assumed  to  be  a  barrier  to  the  Lg  wave  and  the  sampling 
area  is  taken  to  be  bounded  by  that  boudary  ratiier  than  by  the  oceanic 
portion  of  the  ellipse. 

Figure  3.  Stacked  spectral  ratios  from  which  (Q  at  1  hz)  and  the  fre¬ 
quency  dependence  of  Lg  coda  Q  are  obtained.  Those  parameters  are 
determined  by  linear  regression  at  frequencies  between  0.2  and  2.0  hz. 


Figure  4. 
Figure  5. 
Figure  6. 
Figure  7. 
Figure  8. 


Lg  coda  Qo  values  obtained  in  this  study. 

Error  distribution  for  the  image  of  variations. 

Frequency  dependence  of  Lg  coda  Q  at  frequencies  near  1  Hz. 
Error  distribution  for  the  image  of  ti  variations. 

Point  spreading  function  plots  for  four  related  grids. 
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Chapter  4 

Attenuation  of  Surface  Waves  Across 
the  Arabian  Pennisula 

by 

Dogan  Seber  and  Brian  J.  Mitchell 
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ABSTRACT 


Amplitude  spectra  of  fundamental-mode  Rayleigh  and  Love  waves  are  used 
to  determine  shear-wave  Q  (Qp)  models  for  the  upper  crust  of  the  Arabian  pen¬ 
insula.  The  observed  spectra  are  compared  to  theoretical  spectra  computed  for 
earthquakes  with  known  fault-plane  solutions.  Shear-wave  Q  in  the  upper  crust 
can  be  obtained  by  trial-and-error  modeling  to  fit  the  the  slope  of  the  amplitude 
spectra  at  short  periods  (5-30  s),  if  the  source  time  function  is  less  than  about  6  s 
and  if  the  depth  of  the  earthquake  is  known  to  within  about  3  km.  Shear- velocity 
Q  values  obtained  in  this  way  vary  from  about  60  along  the  margin  of  the  Red 
Sea  to  100-150  in  the  central  part  of  the  peninsula,  to  65-80  in  the  eastern  folded 
region.  These  values  are  unusually  low  for  a  stable  region  and  suggest  that 
observed  anelasticity  there  is  affected  by  tectonic  and  epeirogenic  activity  which 
has  occurred  over  the  past  30  m.y.  The  shear-wave  Q  models  derived  from  sur¬ 
face  waves  at  periods  of  5-30  s  closely  predict  reported  values  of  Lg  Q  at  1  Hz; 
shear-wave  Q  therefore  appears  to  be  independent  of  frequency  in  this  region  at 
periods  between  about  1  and  30  s 

INTRODUCTION 

The  Arabian  peninsula  is  a  relatively  small  stable  continental  plate  sur¬ 
rounded  by  two  active  spreading  centers  and  a  zone  of  continent-continent  colli¬ 
sion.  The  crusts  of  most  shields  and  stable  platforms  are  characterized  by  low 
attenuation  of  seismic  waves  (  e.g.  Mitchell,  1973;  Hwang  and  Mitchell,  1987) 
whereas  high  rates  of  attenuation  characterize  spreading  centers  (e.g.  Canas  et 
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al.,  1980;  Chan  et  al.,  1989)  and  zones  of  plate  convergence  (Hwang  and 
Mitchell,  1987;  Correig  et  al.,  1990).  It  is  therefore  interesting  to  study  seismic 
wave  attenuation  in  the  Arabian  peninsula  to  determine  if  it  is  low,  as  it  should 
be  if  it  is  typical  of  that  for  other  shields;  or  if  it  is  high,  as  it  would  be  if  crustal 
anelasticity  there  is  influenced  by  tectonic  process  along  the  borders  of  the  Ara¬ 
bian  plate. 

This  study  presents  determinations  of  seismic  Q  for  the  upper  crust  of  the 
Arabian  peninsula.  Various  paths  across  the  peninsula  allow  us  to  investigate 
the  regional  variation  of  Q  to  learn  whether  or  not  paths  close  to  the  Red  Sea 
spreading  center  differ  in  their  attenuation  characteristics  from  those  further  to 
the  east.  As  a  by-product  of  this  research,  dispersion  of  Rayleigh  and  Love 
waves  is  determined  for  several  paths  across  the  peninsula  and  is  inverted  for 
the  shear  wave  structure  of  the  crust  and  upper  mantle. 

Tectonic  and  Structral  Setting 

I 

The  geodynamic  processes  which  are  active  in  the  Red  Sea  affect  the  present 
geology  and  tectoiucs  of  western  Arabia  (O-Isa  and  Shanti,1989).  Geological  and 
geophysic2il  studies  indicate  that  the  sea  is  still  opening  and  that  new  ocezmic 
crust  is  being  formed  there  (Girdler  and  Styles, 1974;  Girdler  and  Underwood 
1985  ;  Mooney 

et  al.  1985).  An  area  to  the  south  of  the  Peninsula,  in  the  Gulf  of  Aden,  is  also 
a  spreading  center.  Both  the  Red  Sea  and  the  Gulf  of  Aden  spreading  centers 
are  cut  by  many  transform  faults  oriented  in  a  NE-SW  direction  (El-Isa  and 
Shanti  1989  ;  Fairhead  and  Girdler  1970).  The  Arabian  Plate  is  colliding  with  the 
Turkish  plateau  in  the  north,  and  with  the  Persian  plateau  in  the  east  (McKen¬ 
zie,  1972).  Tectonically,  the  Arabian  Peninsula  can  be  divided  into  three  units 
(Knopoff  and  Fouda,  1975).  The  first  unit  is  the  Arabian  shield  which  covers 
almost  one  forth  of  the  entire  peninsula  (Figure  1).  The  second  unit  is  the  stable 
shelf  which  surrounds  the  shield,  and  the  third  unit  contains  the  folded  zones. 
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The  first  unit,  which  is  composed  of  metasediments,  granitic  rocks,  basic  and 
ultrabasic  rocks  and  volcanic  rocks  is  crossed  by  many  faults  that  trend 
northwestward  (Brown,  1970).  The  second  unit,  the  stable  shelf,  is  crossed  by 
major  faults.  The  third  unit  is  the  region  where  most  of  the  folding  occurs  in 
the  peninsula. 

The  United  States  Geological  Survey  and  the  Directorate  General  of  Mineral 
Resources  of  Saudi  Arabia  conducted  a  deep  seismic  refraction  profile  to  obtain  a 
crustal  and  upper  mantle  velocity  model  for  the  shield.  The  results  of  this 
research  show  that  P-wave  velocities  in  the  crust  of  the  Arabian  Shield  are  about 
6.3  km/s  in  the  upper  crust,  about  7.0  km/s  in  the  lower  crust,  and  between  8.0 
and  8.2  km/s  in  the  upper  mantle  (Mooney  et  al.,  1985).  The  average  crustal 
thickness  is  about  40  km  but  thins  to  20  km  in  the  coastal  region  near  the  Red 
Sea.  Niazi  (1968)  obtained  an  average  crustal  and  upper  mantle  velocity  model 
using  phase  and  group  velocity  dispersion  curves  for  the  path  between  stations 
AAE  (Addis  Ababa)  and  SHI  (Shiraz).  His  estimation  for  the  shear  velocities 
was  3.4  km/s  for  the  upper  crust,  3.6  km/s  for  the  lower  crust  and  4.6-4.9  km/s 
for  the  upper  mantle.  Knopoff  and  Fouda  (1975)  measured  the  phase  velocity  of 
Rayleigh  waves  across  the  Arabian  peninsula  using  three  different  profiles. 
They  found  that  there  is  a  low-velocity  channel  for  S-waves  throughout  the 
region  at  a  depth  of  100-140  km.  In  the  charmel  the  velocities  are  about  4.25-4.45 
km/s.  Just  above  this  low  velocity  zone  the  shear  velocities  are  about  4.55  or 
4.65-4.75  km/s.  Ghalib  (1987)  also  determined  velocity  models  from  Rayleigh 
wave  dispersion  curves.  His  results  indicate  shear  velocities  of  about  2-3  km/s 
near  the  surface,  increasing  to  about  5  km/s  at  a  depth  of  100  km. 

Surface  Wave  Attenuation  Studies 

No  measurements  of  attenuation,  from  which  intrinsic  values  of  Q  for  the 
deep  crust  and  upper  mantle  can  be  inferred,  have  been  made  for  the  Arabian 
peninsula.  Ghalib  (1987),  however,  determined  Lg  wave  attenuation  coefficients 
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and  their  frequency  dependence  in  the  peninsula  and  surrounding  areas.  His 
results  show  that  Lg  attenuation  values  in  the  region  at  frequencies  near  1  Hz, 
vary  between  150  and  250  and  the  frequency  dependence  of  Lg  Q  is  about  0.6- 
0.7.  Mokhtar  (1987)  determined  the  attenuation  characteristics  in  the  Arabian 
Shield  at  very  shallow  depths  (  <  .5  km  )  using  short-period  Rayleigh  waves 
obtained  from  a  deep  seismic  refraction  profile  transversing  the  shield.  He 
found  that  at  shallow  depths  values  vary  between  30  and  150. 

Various  methods  have  been  used  to  determine  Q  at  greater  depths  using 
seismic  surface  waves.  These  include  the  two-station  method  where  the  l^-off 
of  spectral  amplitudes  is  measured  directly  at  two  stations  along  a  great  circle 
path  from  an  earthquake  (Tsai  and  Aki,  1969),  radiation  pattern  fitting  methods 
where  observed  spectral  amplitudes  are  compared  with  theoretical  amplitudes 
predicted  by  earthquakes  with  known  fault-plane  solutions  (Tsai  and  Aki,  1969), 
and  methods  where  shapes  of  amplitude  spectra  are  compared  with  predicted 
spectra  (Cheng  and  Mitchell,  1981;  Chan  et  al.,  1989).  For  the  last  of  these 
methods,  Cheng  and  Mitchell  (1981)  utilized  fundamental-mode  and  higher¬ 
mode  Rayleigh  waves,  whereas  Chan  et  al.,  (1989)  utilized  only  fundamental¬ 
mode  waves,  but  were  able  to  use  both  Rayleigh  and  Love  waves  to  advantage. 

Two-station  paths  were  not  available  in  the  present  study,  nor  was  the 
azimuthal  coverage  adequate  for  utilizing  radiation  pattern  methods.  The 
method  remaining  for  us  is  that  of  fitting  spectral  shapes.  Because  higher  modes 
were  poorly  recorded  on  the  seismograms  available  to  us,  we  used  the  method 
of  Chan  et  al.  (1989)  where  spectral  fits  were  sought  for  fundamental-mode  Love 
and  Rayleigh  waves. 

DATA  ACQUISITION  AND  PROCESSING 

Data  for  this  study  consist  of  surface  waves  recorded  on  long-period  instru¬ 
ments  at  three  WWSSN  seismic  stations,  SHI(Shiraz),  TAB(Tabriz),  and  JER 
(Jerusalem).  Five  earthquakes,  all  with  known  fault-plane  solutions,  were 
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employed,  four  of  which  occurred  m  the  gulf  of  Aden,  and  one  in  the  Red  Sea 
(Figure  1).  The  location  and  fault  parameters  for  these  events  are  given  in  Table 
1.  The  seismograms  were  digitized  at  1  $  intervals  and  the  NS  and  EW  com¬ 
ponents  were  rotated  to  obtain  the  transverse  component  of  ground  motion.  To 
separate  noise  and  other  modes  from  the  fundamental  mode,  the  seismograms 
were  phase-match  filtered  using  the  method  developed  by  Herrin  and 
Goforth(1977).  Surface  wave  group  velocities  and  amplitudes  were  obtained  from 
a  multiple  filter  analysis  (Dziewonski  et  al.,  1969)  using  the  vertical  component 
for  Rayleigh  and  the  transverse  component  for  Love  waves. 

Velocity  Models 

Shear  velocity  structure  was  determined  in  two  regions,  western  and 
eastern,  of  the  peninsula.  Average  Love  and  Rayleigh  wave  group  velocity 
dispersion  curves  were  obtained  for  these  two  regions  and  inverted  for  velocity 
structure. 

The  inversion  program  used  in  this  study,  coded  by  D.  Russell  and  R.  B. 
Herrmann,  makes  use  of  partial  derivatives  with  respect  to  layer  parameters. 
Group  and  phase  velocities  are  first  calculated  for  an  initiad  model  and  partiad 
derivatives  of  phase  and  group  velocities  with  respect  to  shear  and  compres- 
sional  velocities  are  obtained  numerically.  Compressional  velocities  are  not 
determined  directly  in  the  inversion  procedure,  but  are  calculated  from  shear 
velocities  by  assuming  a  value  for  Poisson's  ratio. 

An  inversion  procedure  given  by  Lawson  and  Hanson  (1974)  is  used  to 
find  the  unknown  parameters.  If  x  is  the  vector  containing  the  unknown  param¬ 
eters,  A  is  the  matrix  containing  the  partial  derivatives  of  U  with  respect  to  the 
unknown  parameters,  and  b  is  the  difference  between  observational  and  theoret¬ 
ical  values,  then  the  least-squares  solution  to  this  problem  is 


(Ax-b)^  (Ax-b)+a^  x^x=  min. 


(1) 
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where  a  is  an  adjustable  damping  parameter  which  determines  the  tradeoff 
between  resolution  and  standard  deviations  in  the  model  and  can  be  set  to 
achieve  a  desired  degree  of  smoothness.  Matrix  A  can  be  written  in  terms  of 
multiplication  of  three  matrices  by  using  a  singular  decomposition  technique. 

A=U  A  (2) 

U  and  V  matrices  satisfy  the  relation 

UUT=UTU=WT=VTV=I  (3) 

where  I  is  the  identity  matrix.  Then  the  unknown  parameter  vector  can  be  writ¬ 
ten 


x=V{  A2+  I)->  AUT  b,  (4) 

the  variance-covariance  matrix  C  is 

C»V(AH  A2(A2+CT2l)riVT'  (5) 

and  the  resolution  matrix  R  is 

R«V(A2+a2i)-iA^T  (6) 

In  the  inversion  procedure  a  15-layer  model  was  used  where  the  maximum 
depth  of  interest  was  taken  to  be  120  km.  By  using  a  large  number  of  layers  and 
underdetermined  approach  we  avoid  forcing  features  into  a  model  which  might 
be  obtained  using  layers  which  are  thicker  than  those  in  the  Earth.  Smoothness 
can  be  forced  into  the  model  by  selecting  a  suitable  value  for  the  damping  factor 
CT.  The  ratio  of  compressional  velocity  to  shear  velocity  has  been  held  at  V3  , 
which  implies  that  Poisson's  ratio  is  0.25.  Observed  and  predicted  group  velocity 
dispersion  curves,  and  the  corresponding  shear  velocity  structures  for  the 
western  and  eastern  part  of  the  peninsula  axe  given  in  Figures  2  and  3.  As  seen 
from  Figure  2,  standard  deviations  for  Love  and  Rayleigh  wave  group  velocities 
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are  greater  for  the  western  region.  This  is  probably  due  to  the  complex  tectonic 
structure  there.  The  western  peninsula  velocity  model  (Figure  2)  includes  a 
low-velocity  zone  starting  at  a  depth  of  about  70  km.  The  resolving  kernels, 
however,  are  poor  at  this  depth.  Since  the  Red  Sea  is  an  active  spreading 
center,  a  low-velocity  zone  under  that  part  of  the  peninsula  is  not  imexpected. 
The  eastern  region  does  not  include  a  significant  low-velocity  zone  in  the  same 
depth  range.  In  the  eastern  region  shear  velocities  in  the  crust  are  lower  than 
those  of  the  western  region.  Because  of  the  abssence  of  a  low-velocity  zone  in 
the  eastern  model,  they  are,  however,  higher  than  those  int  he  western  region  at 
depths  greater  than  70  km. 


THEORETICAL  SPECTRA 

Theoretical  spectra  were  calculated  using  the  formulation  of  Haskell 
(1963,1964)  as  programmed  by  Wang  (1981).  The  calculations  require  knowledge 
of  the  earthquake  focal  parameters  and  the  velocity-density-Q  model  along  the 
path  of  propagation.  The  focal  parameters  include  the  depth,  strike  and  slip  of 
the  fault  as  well  as  the  source  time  function  and  moment.  The  effects  of  all 
these  parameters  except  moment  on  amplitude  spectra  will  be  examined  in  the 
following  peiragraphs.  Moment  is  neglected  because  it  affects  only  the  level  and 
not  the  slope  of  the  spectra. 

To  study  the  effects  of  different  source  and  fault  parameters  on  the  surface 
wave  amplitude  spectra,  we  performed  a  sensitivity  analysis.  The  spectra  were 
calculated  for  a  source  having  a  depth  of  5  km,  a  dip  of  30“,a  slip  of  30°,  a  strike 
of  30°,  a  four-second  parabolic  source-time  function  and  a  moment  of  5*10^^ 
dyne-cm,  using  the  eastern  peninsula  model.  Spectra  were  computed  for  an 
imaginary  station,  2000  km  away  from  the  source  and  at  an  azimuth  of  0°.  The  Q 
values  were  taken  to  be  150  for  the  first  20  km  and  300  for  the  underlying  half 
space.  Spectra  were  computed  for  cases  in  which  each  parameter  was  allowed 
to  vary  through  a  range  of  values  while  holding  the  others  fixed  at  the  above 
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values. 


The  effect  of  focal  depth 

The  effect  of  the  focal  depth  on  the  spectra  has  been  studied  by  many 
authors  (Tsai  and  Aid,  1970,1971;  Ben-Menahem  and  Sing,1981).  In  this  study, 
the  effect  of  focal  depth  on  the  Rayleigh  and  Love  wave  spectra  is  discussed  for 
focal  depths  of  5  km,  10km,  15  km,  and  20km.  The  amplitude  spectra  for  these 
synthetics  are  given  in  Figure  4.  As  the  depth  of  the  focus  increases,  the  energy 
at  the  shorter  periods  is  reduced  for  both  Rayleigh  and  Love  waves.  The  Love 
wave  spectrum  for  the  longer  periods  is  much  less  sensitive  to  the  focal  depth 
than  it  is  for  shorter  periods.  The  Rayleigh  wave  spectra,  however,  are  sensitive 
to  the  focal  depth  at  both  short  and  long  periods. 

The  effect  of  source  time  function 

The  duration  of  the  source  time  function  is  important  because  it  will  affect 
the  shorter  periods  of  the  Rayleigh  and  Love  wave  spectra.  A  parabolic  source 
time  function  with  different  durations  has  been  used  to  calculate  the  amount  of 
the  effect.  In  Figure  4,  the  Love  and  Rayleigh  wave  spectral  amplitudes  calcu¬ 
lated  for  1  sec,  2  sec,  4  sec  and  6  sec  durations  are  given.  As  the  duration 
increases,  the  effect  of  source  time  function  also  increases  at  shorter  periods. 
The  effect  at  greater  periods  is  insignificant. 

The  effect  of  fault  parameters 

Errors  in  the  determination  of  the  fault  parameters  using  P-wave  first 
motion,  centroid  moment  tensors,  or  surface  wave  amplitude  spectra  are 
thought  to  be  ±10°  for  each  parameter.  To  see  these  effects  on  the  synthetic 
spectra,  spectral  amplitudes  have  been  calculated  for  a  fault  with  dip=30°,  slip= 
30°,  strike  =  30°,  and  ±10°  changes  in  each  parameter. 

Dip 
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The  left  side  of  Figure  5  shows  the  effect  of  variations  in  dip  angle  on  the 
amplitude  spectra  for  Rayleigh  and  Love  waves.  The  effect  of  ±10“  changes  in 
the  dip  angle  is  seen  at  all  period  ranges.  For  Love  waves,  the  spectrum  changes 
in  level  more  than  in  shape.  Rayleigh  wave,  however,  changes  its  shapes  espe¬ 
cially  at  the  longer  periods. 

Strike 

The  center  diagrams  in  Figure  5  shows  the  effect  variation  of  strike  on  the 
amplitude  spectra.  The  Rayleigh  wave  spectrum  varies  for  periods  greater  than 
15  s.  Below  15  s  there  is  almost  no  change  on  the  spectral  shape.  The  Love  wave 
spectra  vary  over  the  entire  period  range  between  5  and  50  s,  with  maximum 
changes  occurring  at  a  period  of  about  20  s. 

Slip 

The  right  side  of  Figure  5  shows  the  effects  of  variations  in  slip  angle  on 
the  spectral  amplitudes.  The  Love  wave  amplitude  spectra  are  affected  more 
than  those  of  Rayleigh  waves,  but  the  change  is  very  small  for  both  cases. 

These  computations  indicate  that  the  spectra  of  both  Love  and  Rayleigh 
waves  are  affected  by  changes  in  the  fault  parameters,  with  most  significant 
changes  occurring  at  the  longer  periods.  The  slopes  of  the  shorter  period  por¬ 
tions  of  the  spectra  are  not  sensitive  to  changes  in  the  fault  parameters.  Since, 
as  we  show  in  the  following  section,  the  effect  of  crustal  Q  values  which  we  are 
interested  in,  is  dominant  at  short  periods,  errors  in  fault  parameters  will  not 
have  adverse  effects  on  our  results. 

The  effect  of  Q 

The  effect  of  Q  on  the  amplitude  spectra  is  investigated  in  this  section.  A 
two-layer  Q  model  is  used  in  which  the  upper  layer  with  a  thickness  of  20  km  , 
overlies  a  half  space.  The  Q  values  used  here  are  assumed  to  be  independent  of 
firequency.  The  effect  of  some  changes  in  the  half  space  Q  values  on  the  spectra 
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are  first  investigated.  The  shear  wave  (Q^)  values  in  the  half  space  were 
adlowed  to  vary  between  150  and  1000  and  the  Q  value  in  the  layer  was  held 
constant  at  150.  The  left  side  of  Figure  6  shows  the  variations  in  the  amplitude 
spectra  of  Rayleigh  and  Love  waves.  The  half  space  Q  values  have  no  significant 
effect  on  the  amplitude  spectra  at  the  shorter  periods.  At  the  longer  periods 
Rayleigh  waves  are  affected  more  than  Love  waves,  but  the  changes  are  rela¬ 
tively  small.  In  the  second  step  values  in  the  layer  were  changed  while  the 
half  space  Q  value  was  held  constant  at  300.  The  right  side  of  Figure  6  shows 
the  crustal  Q  effect  on  the  Rayleigh  and  Love  wave  spectra.  As  seen  from  the 
figure,  Qp  values  for  the  upper  crust  have  a  large  effect  on  the  amplitude  spectra 
at  the  shorter  periods.  The  value  in  the  upper  crust,  in  fact,  plays  the  most 
important  role  in  shaping  the  Rayleigh  and  Love  wave  spectra. 

CALCULATION  OF  Q  VALUES 

We  have  seen  that  the  source-time  function,  the  depth  of  focus  and  Q 
values  in  the  upper  crust  are  most  important  parameters  which  affect  the  shorter 
period  portion  of  the  spectra.  If  therefore,  we  can  constrain  the  depth  of  focus 
and  the  duration  of  the  source-time  function,  we  should  be  able  to  use  the  slope 
of  the  shorter  period  portion  of  the  amplitude  spectra  to  obtain  Q  vcdues  in  the 
upper  crust.  We  use  a  methodology  which  is  similar  to  that  of  Chan  et  al.  (1989) 
in  which  Q  in  the  upper  crust  is  allowed  to  vary  while  all  other  parameters  are 
held  at  reasonable  values.  Chan  et  al.  (1989)  used  spectral  shapes  to  obtain  Q 
values  in  the  Iceland  Plateau.  Since  the  crustal  thickness  is  very  small  near  an 
oceanic  ridge,  they  model  asthenospheric  Qp  values.  Their  procedure  differs 
from  the  present  study  because  they  try  to  match  the  longer  period  portion  of 
the  amplitude  spectra  whereas  we  use  the  shorter  period  portion. 

The  depth  of  each  event  should  be  known  fairly  well  when  using  this  tech¬ 
nique.  The  error  limits  should  not  exceed  ±3  km  to  obtain  a  reliable  Q  values. 
The  earthquakes  that  are  used  in  the  present  study  are  mid-ocean  ridge  type 
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earthquakes,  for  which  the  focal  depths  are  usually  less  than  15  km.  In  Table  1, 
the  source  and  fault  parameters  are  given  for  the  earthquakes  used.  In  the  for¬ 
ward  modeling  we  started  with  an  initial  depth  obtained  from  the  ISC  and  fault 
parameters  obtained  from  centroid  moment  tensor  solutions.  Necessary  changes 
were  made  whenever  the  match  between  observed  and  theoretical  spectra  were 
not  satisfactory  at  the  longer  periods.  Since  the  Rayleigh  wave  amplitude  spec¬ 
tra  are  affected  over  almost  all  periods  with  changing  depth,  it  will  not  be  pos¬ 
sible  to  obtain  a  match  at  the  longer  periods  if  the  depth  is  incorrect,  even 
though  we  have  the  correct  fault-plane  solutions.  Another  constraint  that  we  can 
bring  into  our  calculations  is  to  use  more  than  one  station  for  the  same  earth¬ 
quake.  A  single-station  amplitude  spectrum  may  be  matched  by  more  than  one 
set  of  fault  parameters  if  the  depth  of  the  earthquake  is  changed.  Using  two  or 
more  station  records  will  more  strongly  constrain  the  fault  parameters  and  the 
depth  because  the  same  parameters  must  satisfy  all  amplitude  spectra.  The  use 
of  Love  wave  amplitude  spectra  to  constrain  the  depth  of  focus  will  not  be  use¬ 
ful  because  the  depth  effect  is  only  seen  at  the  shorter  periods,  where  they  may 
trade-off  against  the  effect  of  Q.  The  Love  wave  spectra  may  be  helpful,  how¬ 
ever,  in  constraining  the  fault  parameters  since  they  are  very  sensitive  to  those 
parameters. 

The  source-time  function  is  the  least  constrained  unknown  parameter  in 
our  calculations.  There  is  no  work  that  discusses  the  source  characteristics  of 
earthquakes  in  this  region.  However,  there  are  many  studies  that  discuss  charac¬ 
teristics  of  earthquakes  in  other  regions  (  Tsai  and  Aki,  1970;  Pacheco  and 
Nabelek  1988  ;  Chan  et  cd.  1989).  Tsai  and  Aki  (1970  )  state  that  an  event  of  mag¬ 
nitude  less  than  6  and  with  a  fault  length  about  10  km  will  have  a  rise  time  of 
about  1  s.  Pacheco  and  Nabelek  (1988)  determined  the  source  characteristics  of 
Palm  Springs  earthquakes  having  magnitudes  of  about  6.0  and  estimated  dura¬ 
tions  of  those  earthquakes  to  be  about  4-4.5  s.  Chan  et  al.  (1989)  estimated 
source-time  durations  for  mid-ocean  ridge  earthquakes  with  magnitudes  less 
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than  5.1  to  be  about  4.5  s  or  less.  A  general  conclusion  from  those  studies  is 
that  the  duration  of  the  source-time  function  is  less  than  5  s  for  <m  earthquake 
with  a  magnitude  around  5.0  .  The  earthquakes  that  we  use  in  the  present 
study  have  magiutudes  between  4.9  and  5.3  .  For  this  reason,  it  is  not  likely  that 
the  duration  of  source  functions  will  be  higher  than  5  s  for  the  earthquakes  used 
in  this  study.  Figure  4  indicates  that  source  time  functions  which  might  have  a 
significant  effect  on  the  spectra  have  durations  of  6  s  or  greater.  For  this  reason, 
if  the  source  time  duration  is  less  than  six  seconds,  the  effect  on  the  calculations 
will  be  insigruficant.  Since  six  seconds  is  the  minimum  duration  that  might 
affect  our  calculated  Q  values,  uncertainties  in  source  duration  should  not  affect 
our  results.  If  the  duration  of  the  source  time  function  is  greater  than  six 
seconds,  the  actual  Q  values  will  be  higher  than  what  we  obtain.  For  10  s  dura¬ 
tion  of  the  source  time  function  ,  the  estimated  increase  in  the  Q  values  will  be 
about  20.  A  4  s  parabolic  source  time  function  has  been  assumed  for  all  events. 

In  the  determination  of  Q  models  we  have  not  been  able  to  use  Love  wa^'e 
spectral  amplitudes  for  all  stations  due  to  the  high  noise  levels  on  some 
seismograms.  However,  at  least  one  Love  wave  spectrum  for  each  event  was 
used.  The  use  of  Love  spectra  is  important  because  they  help  to  constrain  the 
fault  parameters  and  also  provide  addition  Q  determinations  for  the  same  event 
and  path  as  that  provided  by  Rayleigh  waves. 

Event  #1  is  located  in  the  Gulf  of  Aden  (Figure  1).  For  this  earthquake  we 
have  records  from  two  stations,  SHI  and  TAB.  Because  of  the  high  noise  level 
we  could  not  obtain  usable  Love  wave  spectra  for  station  TAB.  We  therefore 
have  both  Love  and  Rayleigh  wave  spectra  for  station  SHI  and  Rayleigh  wave 
spectra  only  for  station  TAB.  The  total  of  three  spectra  are  sufficient  to  constrain 
the  depth  of  focus  and  the  fault  parameters.  We  begin  the  forward  mod'^ling 
process  using  initial  values  of  fault  parameters  obtained  from  published  centroid 
moment  tensor  solutions  and  an  initial  value  of  depth  obtained  by  the  ISC.  We 
observe  the  match  between  theoretical  and  observed  spectra  at  the  longer 
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periods  and  change  the  model  parameters  and  the  depth  value  until  we  get  a 
satisfactory  fit  for  the  longer  period  portion  of  the  amplitude  spectra.  Then  we 
introduce  the  Q  values  to  our  calculations.  We  change  the  Q  model  until  we  get 
a  good  match  over  the  short  periods.  For  event  #1  our  results  show  that  Q3 
values  of  about  65  for  Rayleigh  waves  and  about  75  for  Love  waves  give  good 
matches  to  the  spectra  at  station  SHI.  TTie  other  path  gives  a  Qp  value  of  about 
75  at  station  TAB  (Figure  7). 

Event  #2  is  also  located  in  the  Gulf  of  Aden  (Figure  1).  The  available  spec¬ 
tra  for  this  event  are  for  Rayleigh  waves  at  SHI  and  for  both  Love  and  Rayleigh 
waves  at  TAB.  The  models  that  give  a  good  fit  to  the  observed  spectra  have  Q3 
values  between  65  and  70  (Figure  7). 

Event  #3  is  located  very  close  to  event  #2  (Figure  1).  The  available  records 
for  this  event  are  for  Love  and  Rayleigh  waves  at  SHI  and  for  Rayleigh  waves  at 
TAB.  Our  models  for  the  path  to  SHI  gives  Q3  values  about  60  and  85,  respec¬ 
tively,  for  Rayleigh  and  Love  waves.  For  the  path  to  TAB  a  Q3  value  of  60  gives 
the  best  match  to  the  observed  Rayleigh  spectrum  (Figure  7). 

Event  #4  is  located  in  the  Red  Sea.  The  available  records  for  this  event  are 
for  Love  and  Rayleigh  waves  at  SHI  and  for  Rayleigh  waves  at  JER.  These  paths 
cover  different  areas  and  their  azimuths  are  also  different  from  the  other  paths. 
The  results  from  our  modeling  give  Q3  values  of  about  100  and  110,  respectively, 
for  Rayleigh  and  Love  wave  paths  to  SHI  and  about  60  for  the  Rayleigh  wave 
path  to  JER  (Figure  8).  As  seen  from  this  figure,  the  spectrum  for  the  path  to 
JER  has  lost  its  much  higher  frequency  content.  That  lost  of  higher-frequency 
energy  is  caused  by  lower  Q3  values  for  this  path.  This  low  Q  value  may  be 
associated  with  partial  melting  near  the  axis  of  the  Red  Sea  spreading  center 
since  the  path  is  almost  parallel  and  very  close  to  the  ridge  axis. 

Event  #5  was  recorded  only  at  station  JER.  Well-determined  Love  and  Ray¬ 
leigh  wave  spectra  are,  however,  available  for  that  station.  Since  there  is  orUy 
one  station  record,  fault  parameters  and  depth  estimation  for  this  event  are  not 
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well  constrained.  Matching  the  spectra  gives  a  Qp  value  about  150  for  the  upper 
crust  in  the  central  part  of  the  Arabian  peninsula  (Figure  8). 

The  results  obtained  from  forward  modeling  show  that  crustal  values 
across  the  Arabian  Peninsula  vary  between  60  and  150.  In  the  coastal  regions  Qp 
values  are  smaller  than  those  of  inland  regions  in  the  peninsula  (Figure  9). 


Predicted  Lg  Q  values. 

In  order  to  compare  the  results  obtained  in  the  present  study  with  results 
obtained  by  Ghalib  (1987),  synthetic  Lg  waves  were  generated  using  the  western 
velocity  and  Q  models.  Synthetic  Lg  waves  were  generated  for  a  source  at  0.5 
km  depth  and  recordings  have  been  obtained  from  an  imaginary  station  at  dif¬ 
ferent  distances  with  an  azimuth  of  0°.  The  source  parameters  for  that  source  are 
30°  dip,  30°  slip,  and  30°  strike.  By  assuming  that  the  relation 

A(r)=  Ao  e-^  (7) 


where  y  is  the  attenuation  coefficient,  r  is  epicentral  distance,  and  Aq  is  the  ini¬ 
tial  amplitude,  it  is  possible  to  find  Lg  y  by  measuring  the  fall-off  of  Lg  ampli¬ 
tudes  with  distance  (Mitchell  and  Hwang, 1987).  Then  Lg  Q  is  found  using 


where  f  is  frequency  and  U  is  group  velocity. 

In  the  synthetic  Lg  Q  calculations,  U  was  taken  to  be  3.5  km/s.  Calculations 
give  a  predicted  Lg  Q  value  of  about  237  for  this  region,  a  value  which  is  in 
good  agreement  with  the  average  observed  Lg  Q  value  of  230  ±  50  obtained  by 
Ghalib  (1987). 


DISCUSSION  AND  CONCLUSIONS 

Crustal  Qp  values  across  the  Arabian  peninsula  have  been  obtained  using 
surface-wave  amplitude  spectra  for  both  Love  and  Rayleigh  waves.  Theoretical 
amplitude  spectra  calculated  for  a  Q  model  consisting  of  a  single  layer  over  a 
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half-space  gives  a  satisfactory  fit  to  the  observed  data.  The  Qp  values  obtained 
here  include  both  intrinsic  and  scattering  effects.  Since  they  are  relatively  low, 
however,  it  is  likely  that  intrinsic  effects  predominate  (Hwang  and  Mitchell, 
1987). 

Shear  wave  velocity  models  were  also  obtained  from  the  inversion  of 
Rayleigh  and  Love  wave  group  velocity  dispersion  curves.  The  results  from  tiie 
inversion  show  that  shear  velocities  are  about  2.5-3  km/s  in  the  sediments  and 
increase  to  4-4.5  km/s  in  the  deep  crust.  The  results  indicate  a  low-velocity  zone 
in  the  western  part  of  the  peninsula  but  not  in  the  eastern  part.  The  depth  to 
the  top  of  that  low-velocity  layer  is  about  70  km. 

Our  results  pertain  to  average  Qp  values  in  the  upper  crust  of  the  Arabian 
peninsula.  Qp  values  in  the  half-space  are  not  determined  since  this  method  is 
not  sensitive  to  Q  values  in  that  depth  range.  Average  upper  crustal  Qp  values 
are  about  60  in  the  western  coastal  region,  100-150  in  the  central  region,  and 
about  60-85  in  the  eastern  region.  Qp  values  in  all  parts  of  the  Arabian  penin¬ 
sula  are  lower  than  those  in  other  stable  regions  such  as  the  eastern  United 
States  and  Canadian  shelf.  Such  low  values  are  normally  expected  only  in  areas 
having  undergone  relatively  recent  tectonic  activity  (Mitchell,  1975;  Hwang  and 
Mitchell,  1987). 

Bohannon  et  al.  (1989)  described  the  development  of  the  western  portion  of 
Saudia  Arabia  as  beginning  with  volcanism  and  rifting  in  the  Red  Sea  30-32  m.y. 
ago,  followed  by  rotational  block  faulting  and  detachment  faulting  begiiming 
about  25  m.y.  ago,  and  uplift  (2.5  -  4.0  km)  beginning  13.8  m.y.  ago.  The  low  Q 
values  located  in  western  Saudia  Arabia  are  most  likely  due  to  that  recent 
activity.  Relatively  low  values  in  central  Saudia  Arabia  suggest  that  the  activity 
near  the  Red  Sea  had  pervasive  effects  across  much  of  the  peninsula. 

The  low  values  (60-85)  in  the  eastern  part  of  the  study  region  are  difficult  to 
explain  in  terms  of  activity  near  the  Red  Sea  since  it  is  more  than  500  km  from 
the  spreading  axis.  One  possibility  that  slow  deformation,  which  may  be 
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associated  with  past  and  current  tectonic  activity  in  the  Zagros  thrust  zone,  pro¬ 
duces  those  low  values.  Another  possibility  is  that  the  low  Q3  values  there  are 
caused  by  the  thick  blanket  of  sediments  which  covers  the  southeastern  part  of 
the  Arabian  peninsual.  Total  sedimentary  thicknesses  reach  6000  meters  and 
depths  to  the  base  of  Mesozoic  sediment  reach  4500  meters  in  that  region 
(Brown,  1972).  Mitchell  and  Hwang  (1987)  showed  that  young  sediments  could 
cause  a  significant  reduction  in  Q  values  for  surface  waves  and  Lg  in  stable 
regions. 

These  Q3  models  predict  values  of  Qtg  which  are  very  close  to  those  (230 
±50)  obtained  by  Ghalib  (1987).  This  correspondence  indicates  that  Qp  is  not 
only  low  in  this  region,  but  it  is  very  nearly  independent  of  frequency  over  the 
period  range  1-30  s.  If  frequency  dependence  did  occur,  observed  Lg  Q  values  at 
1  Hz  would  be  higher  than  those  predicted  by  the  model  obtained  from  surface 
waves  (Mitchell  1980). 
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Table  1.  Earthquake  soxirce  parameters. 


Event 

Date 

Lat. 

Lon. 

Ori.  time 

Depth 

Mag.  (Ms) 

Dip 

Slip 

Strike 

#1 

02/28/77 

14.78 

54.94 

08:43:55.6 

6  km 

4.9 

70 

-105 

270 

#2 

xnnm 

13.19 

51.01 

23:57:55.0 

4  km 

5.0  Mb 

80 

-85 

280 

#3 

02/11/78 

13.20 

50.98 

12:54:18.0 

9  km 

5.3 

65 

-115 

290 

#4 

16.52 

40.28 

15:00:31.0 

3  km 

5.3 

50 

170 

5 

#5 

01/28/84 

13.96 

51.74 

22:47:50.8 

7  km 

5.3 

70 

150 

150 
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Figure  Captions 


Figure  1.  Tectonic  units  (After  Knopoff  and  Fauda,  1975),  location  of  the  earth¬ 
quakes  and  stations,  and  path  coverage  in  the  Arabian  peninsula. 

Figure  2.  Western  peninsula  dispersion  curves  and  shear  velocity  model.  The 
solid  squares  are  observed  velocities.  The  solid  lines  are  theoretical  dispersion 
curves.  Resolving  kernels  for  selected  depths  are  given  below  the  velocity 
model. 

Figure  3.  Eastern  peninsula  dispersion  curves  and  shear  velocity  model.  The 
solid  squares  are  observed  velocities.  The  solid  lines  are  theoretical  dispersion 
curves.  Resolving  kernels  for  selected  depths  are  given  below  the  velocity 
model. 

Figure  4.  Amplitude  spectra  of  Love  waves  and  Rayleigh  waves  calculated  for 
various  source  depths  (left)  and  source  durations  (right). 

Figure  5.  Amplitude  spectra  of  Love  and  Rayleigh  waves  calculated  for  various 
dip  values  (left),  strike  values  (center),  and  slip  values  (right). 

Figure  6.  Amplitude  spectra  of  Love  and  Rayleigh  waves  calculated  for  various 
values  of  at  depths  greater  than  20  km  (left),  and  at  depths  less  than  20  km 
(right). 

Figure  7.  Theoretical  and  observed  amplitude  spectra  at  SHI  and  TAB  for  event  1 
(top),  event  2  (center),  and  event  3  (bottom).  The  solid  circles  are  observed  spec¬ 
tral  values  and  the  solid  lines  are  the  theoretical  amplitude  spectra  calculated  for 
a  single  layer  over  a  half  space  model.  The  broken  lines  are  the  theoretical 
amplitude  spectra  calculated  for  upper  crustal  values  which  are  greater  than. 


115 


and  less  than,  the  derived  value  by  20. 

Figure  8.  Theoretical  and  observed  amplitude  spectra  at  stations  SHI  and  JER  for 
event  4  (upper)  and  event  5  (lower).  See  the  caption  of  Figure  7  for  explanations 
of  the  symbols. 

Figure  9,  values  in  the  upper  crust  of  the  Arabian  peninsula. 
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Chapter  5 


Learned  Estimation  of  Focal  Depths  with 
Multilayered  Netiral  Networks:  Application 
to  Continental  and  Oceanic  Earthquakes 

by 


Gagan  B.  Patnaik,  Dogan  Seber  and  Brian  J.  Mitchell 


ABSTRACT 

A  novel  technique  based  upon  neural  network  learning  has  been  applied  to  the 
estimation  of  earthquake  focal  depths  from  observed  amplitude  spectral  ratios  of  surface 
waves.  The  network  learning  exploits  mapping  relationships  between  the  spectral  ratios 
of  fundamental-mode  synthetic  Love  to  Rayleigh  wave  spectra  and  corresponding 
source  depths.  On  a  partially  trained  network  with  a  continental  velocity  model,  4608 
noise-free  synthetic  spectral  ratios  corresponding  to  various  source  orientations  were 
consistently  classified  into  8  preassigned  focal  depth  groups  with  over  80%  accuracy 
when  the  networks  were  simulated  with  varying  numbers  of  hidden  layers  and  varying 
numbers  of  processing  units  in  each  layer.  With  this  satisfactory  performance  of  trained 
neural  networks  on  synthetic  spectra,  well-calibrated  events  from  eastern  North  Amer¬ 
ica  were  processed  and  the  corresponding  focal  depths  were  estimated.  The  events  stu¬ 
died  in  this  paper  were  the  New  Brunswick  earthquake  of  09  Jan  1982  and  three  of  its 
aftershocks.  The  data  for  all  events  of  the  New  Brunswick  sequence  were  retrieved 
from  GDSN  stations  within  continental  North  America  (ALQ,  ANMO,  JAS,  LON). 
Phase-matched  filtered  fundamental-mode  surface  wave  signals  were  extracted  for  the 
construction  of  spectral  ratios.  Focal  depths  obtained  in  all  cases  were  consistent  with 
those  reported  by  others. 

A  second  neural  network  with  a  different  architecture  was  trained  using  spectral 
ratios  computed  from  an  oceanic  velocity  model  derived  from  the  inversion  of  surface 
wave  dispersion  across  the  southeast  Indian  Ocean  for  focal  depth  estimation  of  oceanic 
events.  The  network,  after  partial  training  with  3240  synthetic  spectral  ratios,  classified 
them  into  5  preassigned  focal  depth  groups  with  over  90%  accuracy.  Two  events  which 
occurred  on  the  southern  Indian  Ocean  were  chosen  for  testing.  The  depth  of  one  of 
the  events  agreed  within  ±2.5  km  of  the  reported  centroid  depth  as  well  as  the  ISC 
reported  depth  from  body  wave  phases  (pP  -  P).  The  estimated  depth  of  the  other 
event  also  agreed  with  the  depths  assumed  for  that  region  by  the  NEIS.  Since  the  syn¬ 
thetic  spectral  ratios  were  noise-free,  a  noise  tolerance  test  was  also  performed  using  a 
subset  (300)  of  the  spectral  ratios  generated  with  the  oceanic  model.  Upon  training  and 
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testing  at  various  noise  levels  it  was  observed  that  even  with  40%  noise  added  to  the 
spectral  ratios  (a  S/N  degradation  of  8dB),  the  network  was  still  able  to  classify  the  spec¬ 
tral  ratios  into  3  preassigned  depth  groups  at  nearly  the  70%  level.  The  noise  tolerance 
test  implied  that  it  might  even  be  possible  to  use  observed  spectral  ratios  with  low 
signal-to-noise  for  depth  estimation.  Further  testing  with  many  observations  will  be 
necessary  to  determine  whether  or  not  this  will  be  possible. 


Introduction 

Recent  developments  in  the  field  of  neural  computing  have  provided  efficient  non¬ 
linear  transfer  function  approximation  mechanisms  {e.g. ,  Rumelhart,  et  al.  1986;  Qark, 
1988;  Moody  and  Darken,  1989)  that  can  be  applied  to  many  seismological  problems. 
Neural  computing  may  also  be  viewed  as  a  decision-making  process  where  a  neural  net¬ 
work  is  trained  to  identify  the  class  which  may  best  represent  an  input  pattern.  This,  in 
turn,  requires  that  a  network  has  already  been  adapted  to  recognizing  such  classes,  by 
means  of  learning  earlier  examples.  The  learning  can  essentially  be  viewed  as  a  non¬ 
linear  system  modeling  process.  It  so  happens  that  these  networks,  in  some  restricted 
sense,  are  capable  of  performing  certain  low-level  tasks  such  as,  combining  numerical 
information  and  pattern  recognition.  Solutions  to  the  problem  of  the  present  study, 
focal  depth  estimation,  are  sought  by  exploiting  this  capability  of  neural  networks. 

The  depth  dependence  of  surface  wave  amplitude  spectra  has  been  well  studied, 
(e.g. ,  Toksoz  et  al. ,  1964;  Harkrider,  1964,1970;  Tsai  and  Aki,  1970,1971;  Masse  et  al. , 
1973;  Herrmarm,  1979;  Yan  and  Alexander,  1990).  The  underlying  principle  can  be 
described  as  follows.  There  is  a  strong  dependence  of  observed  surface  Rayleigh  wave 
amplitude  spectra  on  the  depth  of  the  source  for  many  focal  mechanisms.  The  spectral 
holes  of  fundamental-mode  Rayleigh  waves  shift  towards  lower  frequencies  as  the  focal 
depth  increases.  In  addition,  the  high-frequency  spectrum  of  Love  and  Rayleigh  waves 
attenuate  as  the  depth  increases  and  the  slope  of  the  spectrum  becomes  steeper  with 
decreasing  cut-off  frequency.  There  is  also  a  strong  azimuthal  dependence  of  the  spec¬ 
trum  for  a  given  source  mechanism.  These  facts  have  been  used  in  the  past  to  estimate 
the  source  parameters  of  shallow  events. 

The  aim  of  this  study  is  to  estimate  focal  depths  from  observed  surface  wave 
amplitude  spectral  ratios  by  processing  those  ratios  through  a  neural  network  that  has 
been  a  priori  adapted  to  learn  the  nonlinear  functional  relationship  that  may  exist  among 
the  spectral  ratios  (input),  source-crustal  response  and  source  mechanisms  (system), 
and  the  focal  depths  (output).  By  using  Love  to  Rayleigh  wave  spectral  ratios,  we  effec¬ 
tively  cancel  out  attenuation  factors  (geometrical  spreading  and  anelasticity).  Incor¬ 
porating  a  well-constrained  source-crustal  model  response  and  combinations  of  source 
mechanisms  into  the  neural  network  learning  process,  variations  in  the  spectral  ratios 
can  be  attributed  to  changes  in  source  depths. 
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In  the  following  sections  we  present  brief  introductions  to  neural  computing,  depth 
estimation  from  surface  wave  amplitude  spectral  ratios,  and  preprocessing  methods 
applied  to  the  observed  data.  We  then  apply  the  present  technique  to  four  events  of 
New  Brunswick  sequence  emd  two  events  in  the  Indian  Ocean  with  reported  focal 
depths. 


Neural  Computing 

Detailed  theoretical  developments  of  neural  networks  have  been  presented  by 
many  researchers,  primarily  in  the  fields  of  electrical  engineering,  physics,  biology  and 
computer  sciences  (e.g. ,  Hopfield,  1982;  Hopfield  and  Tank,  1985;  Rumelhart  and 
McQelland,  1986;  Sejnowski  et  al. ,  1986;  Takeda  and  Goodman,  1986;  Lippmaim,  1987; 
Peterson  and  Anderson,  1987).  A  neural  network  is  shown  in  Figure  1.  It  consists  of 
<N>  layers,  each  containing  S[i]  noidinear  processors  (neurons),  i  =  l,..,n.  Successive 
layers  are  fully  interconnected  in  the  forward  direction  (feed-forward).  This  architecture 
is  computationally  complete  (Lippmann,  1987)  in  the  sense  that  any  arbitrary  classifica¬ 
tion  boundary  c>*n  ’«  approximated  with  this  topology.  Each  processor  acts  as  a  thres¬ 
holding  unit,  tun  ing  "on"  or  "off"  depending  on  whether  or  not  the  sum  of  its  inputs 
exceeds  a  th  ^shold  value.  This  output  'f  a  processor  becomes  an  input  to  the  other 
processors  after  being  multiplied  by  a  set  of  interconnection  weights.  Determination  of 
the  appropriate  interconnection  weights  (strengths)  among  all  the  processors  such  that 
the  dynamic  system  as  a  whole  (the  network)  produces  a  desirable  output  constitutes 
the  network  learning .  In  general,  when  a  network  has  learned ,  an  energy  functional  has 
been  minimized,  and  the  matrix  representing  the  interconnection  topology  yields  a 
trained  neural  network.  This  can  be  achieved  using  a  numerical  iterative  procedure,  by 
solving  a  system  of  stiff  ordinary  differential  equations  representing  the  network 
dynamics. 

In  our  case  a  network  learns  the  given  spectral  ratio  inputs  and  the  corresponding 
focd  depth  outputs  in  a  best  match  sense  and  the  heart  of  this  learning  is  an  iterative 
nonlinear  optimization  process  where  a  large  number  of  weak  constraints  are  simultane¬ 
ously  satisfied.  The  underlying  nonlinear  system  which  is  modeled  is  the  layered- 
medium  response  to  the  matching  of  given  spectral  ratios  and  corresponding  focal 
depths.  There  are  several  ways  that  a  neural  network  can  be  trained.  The  following 
paragraphs  describe  two  of  those  which  are  used  in  this  study. 

The  Back-Error-Propagation  Learning  Algorithm 

This  is  a  recently  developed,  effective  learning  algorithm  for  the  multi-layer  neural 
network  (multilayer  perceptron)  shown  in  Figure  1.  Following  Rumelhart  and  McClel¬ 
land  (1986),  and  Khotanazad  and  Lu  (1990),  learning  can  be  accomplished  by  the  optimi¬ 
zation  of  a  criterion  function,  where  constraint  satisfaction  can  be  achieved  by 
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estimating  the  discrepancy  between  the  desired  and  actual  output  values,  feeding  back 
an  error  signal  layer  by  layer  to  the  inputs,  and  then  adjusting  the  interconnection 
weights  in  such  ways  as  to  modify  them  in  proportion  to  their  contribution  to  the  total 
mean-square-error.  This  is  computationally  the  fastest  deterministic  algorithm;  however, 
it  is  also  the  most  susceptible  to  local  minima  problems.  Following  Khotanazad  and  Lu 
(1990),  this  learning  procedure  can  be  briefly  described  as  follows.  The  connection 
weights  Wjj  are  initialized  to  small  random  values.  A  feature  vector  (spectral  ratio)  is  fed 
forward  through  the  multilayer  neiual  network.  Then  the  output  of  processor  (neuron)  j 
which  is  a  nonlinear  sigmoidal  over  its  total  input  becomes, 

S,  ■  sig(  ) - (1) 

'  1  +  e  ' 


In  (1),  Sj  is  the  output  of  the  i*  neuron  in  a  previous  layer.  The  function  sig(x)  is  of  the 
form  of  tanh  or  of  a  modified  Heaviside  function.  Wjj  and  0j  are  interconriection  weights 
and  thresholds,  respectively.  The  squashing  function  sig(x)  with  a  sigmoidal  shape  con¬ 
trols  the  nonlinear  behavior  of  each  processor.  Next,  an  error  term  Sj  for  a  processor  at 
the  output  layer  is  computed.  That  term  is  given  by 

8,  =  (Oj-Sj)S,(l-S,)  (2) 


where  Oj  is  the  desired  output  controlled  by  an  apriori  teaching  signal,  and  Sj  is  a 
neural  output.  Similarly,  using  the  chain-rule  of  differentiation,  the  error  term  for  a  hid¬ 
den  layer  processor  is  computed, 

8j  =  Sj(l-Sj)25kWjk  (3) 

k 


where  k  runs  over  all  processors  in  the  layer  above  processor  j.  Finally,  the  connections 
are  adjusted  by  using. 


Wjj'*''  =  Wjj""'"’  +  -nSjSj  +  a(Wi|“"'"‘-  Wif*^). 


(4) 


The  parameter  ii  is  a  learning  rate  similar  to  step  size  in  gradient  search  algo¬ 
rithms,  and  a  is  a  constant  that  reflects  the  effect  of  past  weight  changes.  The  latter  is 
somewhat  similar  to  momentum,  a  parameter  that  effectively  filters  out  high-frequency 
variations  of  an  error  surface  and  helps  avoid  local  minima.  Thus,  in  (4)  not  only  the 
gradient,  but  also  the  previous  change  in  interconnections  has  been  used  to  adjust  the 
connection  vectors.  This  is  analogous  to  a  second-order  dynamic  system  that  has  better 
convergence  properties  than  a  first-order  system  (gradient  descent  without  a  term  in 
(4)). 

All  training  input  vectors  are  presented  either  cyclically  or  randomly  to  the  neural 
network,  errors  are  propagated  backwards  and  the  cycles  are  continued  until  the  con¬ 
nection  matrix  stabilizes.  The  nonlinear  squashing  function  used  in  each  processor  may 
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be  replaced  by  other  simpler  or\es  (e.g.,  linear  or  piece -wise  linear  functions),  if  the 
feature  space  is  simple.  This  would  effectively  speed  up  the  learning  process.  In  our 
simulations,  however,  aU  the  processors  have  a  sigmoidal  threshold. 

Statistical  Mechanics  Analogy 

In  neural  network  learning  processes  constraint  satisfaction  is  achieved  through 
large  scede  nonlinear  optimization  such  that,  either  a  functional  describing  the  goodness 
of  fit  is  mckximized,  or  equivalently,  a  functional  describing  the  energy  of  a  dynamical 
system  (state  of  the  network)  is  minimized.  The  other  learning  algorithm.  Mean  Field 
Theory,  used  in  this  study,  utilizes  this  energy  minimization  criteria  from  a  statistical 
physics  perspective.  The  main  advantage  of  this  over  other  deterministic  procedures  is 
its  ability  to  overcome  the  problem  of  local  minima .  In  a  physical  system  the  tendency  is 
always  to  evolve  from  highly  energetic  states  to  states  of  minimal  energy.  From  this 
perspective  the  neural  network  dynamics  may  completely  be  described  as  thermo¬ 
dynamical  systems  that  seek  stable  minimum  energy  states,  or  attractor  states  (Hopfield, 
1982,  1985).  Many  recent  investigations  {e.g. ,  Moussouris,  1974;  Choi  and  Huberman, 
1983;  Geman  and  Geman,  1984;  Sejnowski,  et  al.  1986;  Sompolinsky,  1988;  and  Clark, 
1988),  report  exhaustive  analytical  descriptions  of  neural  networks  through 
attractor  dynamics  and  study  their  behavior  as  deterministic  chaotic.  Most  of  these  stu¬ 
dies  are  of  theoretical  importance  and  we  will  restrict  ourselves  to  only  application- 
oriented  analyses. 

Since  constraint  satisfaction  is  an  energy  minimization  problem,  an  analogy  may  be 
drawn  to  metal  annealing  processes,  which  involve  heating  the  material  and  then  cool¬ 
ing  it  very  slowly,  so  that  the  temperature  settles  gradually  and  the  material  at  molecu¬ 
lar  level  coalesces,  giving  rise  to  a  global  energy  minimum.  Hence,  the  idea  of 
simulated  annealing  of  a  dynamical  system  to  bring  it  to  a  stable  state  which,  in  turn, 
yields  the  appropriate  globed  solution.  In  practice,  it  is  done  through  a  probabilistic 
update  rule, 

where  T  is  a  global  temperature  peirameter  and  the  term  netj  corresponds  to  (1).  In  this 
case,  the  processors  may  become  continuous,  depending  upon  the  probability  distribu¬ 
tion  in  (5)  and  and  the  system  as  a  whole  becomes  stochastic. 

The  Mean-Field-Theory  Learning  Algorithm 

Mean  Field  Theory  (MFT)  approximation  is  a  well  known  technique  in  spin- 
systems  {e.g.,  Glauber,  1963).  Peterson  and  Anderson  (1987)  have  improved  an  earlier 
non-deterministic  learning  algorithm  (the  Boltzmann  Machine)  using  this  approximation 
such  that  the  solutions  to  a  set  of  deterministic  equations  could  be  used  to  describe  the 
dynamic  state  of  a  network.  The  application  of  MFT  to  neural  network  learning  results 
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in  a  statistical  mechanical  system  that  uses  a  simulated  annealing  relaxation  technique 
to  search  for  the  best  possible  connection  weights  (globally  minimum  energy)  over  a 
large  solution  space. 

An  energy  functional  may  be  written  as, 

1  N+h 

E(S)  =  -i2W,S,S|  (6) 

^  i.i-0 

where  S  =  (Sj,  .  .  .  ,  Sj,  ■  •  •  ,5^)  describes  a  state  of  the  network,  Wjj  are  interconnection 
weights,  h  is  the  number  of  hidden  processors,  and  Sj  j  are  binary  valued  processors 
(nevurons).  Differentiating  (6)  with  respect  to  Sj  and  noting  the  symmetric  coupling, 
Wjj  =  Wji,  the  system  can  be  made  to  relax  to  a  local  energy  minima  with  a  step  function 
updating  rule. 

Si  =  +1  ;  SWijSj  >0  (7) 

i 

=  -1  ;  XWijS,  <0. 

i 

The  simulated  annealing  technique  is  applied,  in  practice,  to  avoid  any  local  minima 
and  the  relaxation  is  done  according  to  the  Boltzmann  distribution, 

P(S)  a  e-^S)/T  (g) 

where  P(S)  is  the  probability  that  the  state  S  will  occur  given  the  temperature  T. 

Peterson  and  Anderson  (1987)  derived  an  average  of  a  network  state  function  in  (8)  as, 

<F{S)>  =  4-  2  P(S)  (9) 

^  s 


where  Z  is  the  so-called  partition  function, 

(10) 

S 

The  summation  is  taken  over  all  the  processor  configurations,  therefore,  the  sum  over 
the  binary  processor,  S=±l  in  (7),  can  be  replaced  by  an  integral  over  continuous  mean 
field  variables.  With  further  reductions  the  solutions  to  the  neural  network  represented 
in  (6)  may  be  written  in  their  final  form  as 

Sj  =  tanh(XWijSi/T).  (11) 

i 

Note  that  the  unit  step  function  approximation  used  in  (7)  is  now  replaced  by  a  sig¬ 
moidal  updating  rule.  The  algorithms  based  on  (11)  are  deterministic  and  can  be  solved 
iteratively  by 
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(12) 


:  =  tanh(2WijSi®'‘*/T). 

1 

In  brief,  each  spectral  ratio  (feature  vector)  is  passed  through  two  phases  (clamped 
and  free  running)  each  with  an  annealing  schedule.  Statistics  are  collected  at  the  end  of 
each  phase  and  the  interconnection  weights  are  updated.  The  steps  are  repeated  until 
no  more  changes  in  Wjj  take  place. 

Both  learning  algorithms,  the  back-error-propagation  and  the  mean-field-theory, 
have  been  used  in  the  present  study. 


Spectral  Ratios 

Following  Haskell  (1964)  and  Harkrider  (1964,  1970),  the  spectral  far-field  surface 
displacements  for  Love  and  Rayleigh  waves  due  to  a  double-couple  source  in  a  layered 
earth  may  be  written  as, 

U  =  S  k®  X(0/h)  E  (13) 

where  S  is  the  spectral  source  function,  k  is  the  angular  wave  number  (  k=ayfc  ),  X  is  the 
source  radiation  pattern,  h  is  the  source  depth,  0  is  the  azimuth  from  source  to  station, 
and  r  is  the  epicentral  distance.  (13)  may  be  separated  into  Love  and  Rayleigh  com¬ 
ponent  displacements  as 

Ul  =  S  kL  Xl(0.E)  El  Dl  (kL,r)  (14) 

and 

Ur  =  S  kR  e-^"^  XR(0,h)  Er  Dr  (kR,r).  (15) 

where  Dl  r  (k,r)  represents  attenuation  effects.  This  follows  from  the  fact  that  in  (13), 
m=l  for  a  double-couple  source.  The  quantities  El  and  Er  can  be  determined  from  the 
amplitude  responses,  Al,r/  the  wave  numbers,  kL,R,  and  the  Rayleigh  wave  ellipticity, 
€o,  using  El  =  Al  ^nd  Er  =  Cg  Ar  kR^'*  .  Dl,r  (k,r)  can  be  written  as 

-•rrfr 

Dl,r  (k,r)  -  g  Qu.Uu, 

which  includes  attenuation  due  to  both  geometrical  spreading  and  anelasticity.  Ql,r  in 
this  equation  represents  the  quality  factor  for  either  Love  and  Rayleigh  waves,  and  Ul,r 
is  the  corresponding  wave  group  velocity.  Since  we  are  interested  in  spectral  ratios  of 
Love  to  Rayleigh  waves,  with  the  assumption  that  the  attenuation  coefficients  for  both 
these  waves  are  approximately  the  same  for  the  frequency  range  of  observable  continen¬ 
tal  dispersion  (e.g.,  Tsai  and  Aki,  1970;  Yan  and  Alexander,  1990),  the  ratio 
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DL(kL/r)/DR(kR,r)  is  very  dose  to  1,  Le.,  the  terms  involving  attenuation  and  geometri¬ 
cal  spreading  cancel  out.  The  quantity  in  (13)  that  remains  to  be  defined  is  the  complex 
radiation  pattern  function  given  by 

X(6,h)  =  do  +  i(disin0  +  d2Cos0)  +  d3sin20  +  d4cos20.  (16) 


The  coeffidents  dj,  as  given  in  Harkrider  (1970),  may  be  expressed  in  terms  of 
Thomson-Haskell  displacement-stress  vector  elements  (Haskell,  1953).  func¬ 

tions  of  source  depth  (h),  dip  (8)  and  slip  (X.)  of  the  fault  plane,  and  the  azimuth  (0). 

The  amplitude  responses,  Al  and  Ar  ,  and  hence  the  quantities  El  and  Er  ,  are 
functions  of  frequency,  mode  order  and  the  elastic  properties  of  the  layered  medium 
and  are  independent  of  source  geometry  and  depth.  The  term  S  in  (13)  is  the  source- 
time-space  function  (Harkrider,  1970)  and  will  depend  on  fault  displacement,  fault 
dimensions,  and  rupture  velocity.  For  smaller  sources,  with  magnitudes  less  than  6.0, 
the  assumption  of  non-propagating  point  source  may  be  valid.  Therefore,  ignoring  the 
source  finiteness,  S  becomes  the  same  for  both  Love  and  Rayleigh  waves  and  will  cancel 
out  in  the  ratios  for  periods  greater  than  10  seconds.  Taking  the  ratio  between  (14)  and 
(15)  we  get. 


Ul 

Ur 

XRfO/I'A/S)  Er 

The  first  and  last  quantities  on  the  right-hand-side  of  (17)  are  model  dependent 
whereas,  the  complex  radiation  pattern  terms  are  dependent  on  source-crustal  model, 
the  source  geometry,  and  the  source  depth.  This  formulation  has  been  used  by  many 
for  source  mechanism  and  source  depth  estimation  {e.g.  Harkrider,  1970;  Masse  et  al. , 
1973;  Herrmann,  1979;  Van  and  Alexander,  1990).  In  the  present  study  we  also  use  (17) 
to  construct  the  fundamental-mode  surface  wave  spectral  ratios. 


Data  Processing 

It  is  weU  known  that  observed  surface  wave  modes  are  frequently  contaminated  by 
both  random  and  coherent  noise  due  to  multipathing  and  higher-mode  interference. 
Therefore,  a  direct  Fourier  amplitude  spectrum  of  the  observed  signal  is  not  a  true 
representation  of  the  fundamental-mode  surface  wave  spectrum.  There  are  many  tech¬ 
niques  available  for  the  isolation  of  the  pure  fundamental-mode  (or  any  higher  mode) 
spectrum.  The  phase-matched  filtering  method  has  been  used  by  many  researchers 
{e.g. ,  Herrin  and  Goforth,  1977;  Goforth  and  Herrin,  1979;  Russell,  1987)  for  surface 
wave  studies.  This  method  essentially  involves  compression  of  energy  in  the  time 
domain  so  that  a  zero-phase  signal  with  energy  concentrated  about  zero-lag  results.  The 
resulting  signal  with  no  phase  is  equivalent  to  a  time-domain  pseudo-autocorrelation 
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function.  Then  by  time-windowing  the  pseudo-autocorrelation  function,  il;(t),  with  a 
symmetric,  zero-phase  window,  w(t),  and  taking  a  Fourier  transform,  the  desired  signal 
spectrum  is  isolated  from  other  modes  and  noises.  The  amplitude  spectrum  of  the 
pseudo-autocorrelation  function  may  be  written  as, 

00 

Ae***"  =  /  4»(t)  w(t)  e-‘“‘dt  (18) 

.jao 


where  8kx  =  (fe  -  k)x  represents  r*  residual  phase  error  term.  This  error  will  depend  on 
the  initial  estimate  of  the  wave  number  k  which  can  be  iteratively  unproved  as. 


To  estimate  an  initial  wave  number,  Herrin  and  Goforth  (1977)  used  the  multiple- 
filter  technique.  Many  recent  studies  involving  surface  wave  dispersion  describe  this 
technique  in  detail  (e.g.,  Dziewonski  et  al.,  1969;  Herrmann,  1973;  Bache  et  al. ,  1978; 
and  Yan  and  Alexander,  1990).  Following  Herrmann  (1973),  a  narrow  bandpassed  fil¬ 
tered  seismogram  may  be  constructed  as 


g(t/X)  = 


2ir 


V2 

JL 

a 

\  > 

2A„((o^,,x)  exp  i  [o)„t  -  ko,„x]  exp 

■4“  {'“Uo,*) 

(19) 


where  Ajn(a)j,,x)  is  the  amplitude  spectrum  of  the  m***  mode  at  frequency  and  the 
narrow-band  Gaussian  filter  used  to  construct  (19)  has  the  form. 


H(a))= 


exp(-a(i)Vo)^) 

0 


I  0)  I  <(i)^ 

I  Cl)  I  >(i)j.  , 


The  cut-off  frequency  (coj  and  the  filter  bandwidth  parameter  (a)  have  been  chosen 
in  this  study  to  be  co^  =  cjo/4,  and  a  =  16Tr,  respectively.  Note  that  in  (19)  Ug  n,  is  the 
group  velocity  of  the  chosen  mode,  m,  at  frequency  o),,  which  can  be  written  as, 
U<,m  =  x/tg.  Since  the  epicentral  distance  x  is  known,  from  the  times  corresponding  to 
the  maxima  (group  delay  time)  of  the  modulus  of  (19),  a  group  delay  curve  can  be 
obtained.  From  these  group  delay  times  of  the  signal  at  various  frequencies,  one  can 
construct  a  group  velcxrity  dispersion  curve.  The  application  of  the  multiple-filter  tecdi- 
nique  is  effective  when  the  wave  trains  belonging  to  two  or  more  modes  are  superim¬ 
posed  in  time,  so  long  as  these  superimposed  modes  occur  at  different  frequencies. 
From  the  group  delay  times  of  the  fundamental  mode  signal  at  various  frequencies,  one 
can  construct  a  phase  spectrum  given  by 

<}>=/  tg  do).  (20) 

0 
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This  equation  is  used  to  compute  the  initiad  estimate  of  the  wave  number,  k, 
(<j)  =  kx)  for  the  phase-matched  filter.  In  practice,  a  family  of  wave  numbers  will  result 
corresponding  to  multiples  of  2it,  This  can  be  resolved  with  a  knowledge  of  one  stable 
phase  velocity  estimate  at  a  given  period  for  the  region  under  consideration. 

The  time-domain  windowing  (18)  of  the  pseudo-autocorrelation  function  results  in 
biased  estimates  in  frequency  domain.  Detailed  considerations  for  reducing  the  bias  and 
a  suitable  choice  of  window  functions  used  in  the  present  study  are  given  in  Jenkins 
and  Watts  (1968)  and  RusseU  et  al.  (1988).  A  Parzen  window  is  used  for  processing  the 
phase  spectrum  while  a  cosine-tapered  window  is  used  to  extract  the  amplitude  spec¬ 
trum.  The  total  bias  error  due  to  windowing  can  safely  be  under  5%.  Since  we  are 
interested  in  normalized  spectral  ratios,  such  bias  may  not  be  significant.  In  this  study, 
both  of  these  filtering  techniques  were  utilized  to  extract  the  fundamental-mode  ampli¬ 
tude  spectra  from  each  observed  long-period  seismogram. 


Simulations  and  Observations 

A  schematic  drawing  (Figure  2)  shows  the  overall  methodology  adopted  in  this 
study.  First,  a  suitable  velocity  model  is  chosen  for  the  region  under  study.  Synthetic 
spectral  ratios  of  Love  to  Rayleigh  waves  are  computed  and  a  carefully  designed  neural 
network  is  then  trained  with  these  synthetic  data.  Necessary  tuning  of  the  network  is 
done  during  the  training  and  testing  phase.  If  the  performance  becomes  satisfactory 
upon  testing  with  the  training  data  set,  synthetic  spectral  ratios  at  intermediate  depths 
(other  than  those  chosen  for  training)  are  generated  for  testing.  Next,  observed  surface 
wave  spectral  ratios  are  constructed  from  as  many  recorded  seismograms  as  possible, 
which  when  processed  through  the  above  trained  network  provides  a  focal  depth  esti¬ 
mation  for  that  event.  However,  if  for  a  given  region  well  distributed  observed  events 
are  available  with  pre-estimated  focal  depths,  then  a  network  may  be  constructed  with 
the  corresponding  observed  spectral  ratios.  Such  a  network  should  be  ideally  desirable, 
the  advantage  being  that  inherent  noise,  characteristic  to  the  area  under  consideration, 
would  also  be  incorporated  in  the  learning  process.  Since  such  situations  (enough  obser¬ 
vations  for  an  appropriate  mapping)  are  rarely  true,  our  alternative,  training  a  network 
with  synthetic  data  derived  from  well  constrained  velocity  models  representing  an  area 
may  be  useful  for  focal  depth  estimations. 

Synthetic  Spectral  Ratios 

Well  constrained  shear-velocity  models  for  Eastern  United  States  have  been  derived 
by  many  investigators  {e.g. ,  Mitchell  and  Herrmann,  1979;  Taylor  and  Toksoz,  1982).  In 
this  study,  the  Taylor  and  Toksoz  (1982)  velocity  model  (Figure  3),  derived  from  the 
inversion  of  surface  wave  group  and  phase  velocities,  has  been  used.  This  model  has 
also  been  used  by  Yan  and  Alexander  (1990)  for  the  estimation  of  focal  depth  and  fault 
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orientations  of  the  New  Brunswick  sequence.  Since  one  of  our  aims  is  to  compare  the 
results  obtained  from  our  technique  for  the  aftershock  depths  with  those  of  Yan  and 
Alexander  (1990),  it  is  essential  to  use  the  same  model.  Those  authors  showed  that  the 
velocity  model  explained  the  observed  average  dispersion  along  the  paths  from  the  New 
Brunswick  source  region  to  GDSN  stations  used  in  their  study,  very  well. 

Using  the  normal  mode  summation  technique  (Wang,  1981),  fundamental-mode 
Love  and  Rayleigh  wave  spectra  were  generated  and  the  corresponding  ratios  obtained. 
Since  the  neural  network  training  process  is  very  time  consuming,  we  restricted  the 
source  parameters  as  follows:  depths  were  restricted  to  3  km  and  17  km  with  an  interval 
of  2  km,  dip  ranged  from  10®  to  80®  with  an  interval  of  10®,  slip  ranged  from  10®  to  150®, 
with  an  interval  of  10®  and  strike  ranged  from  175®  to  225°  with  an  interval  of  10®.  The 
fault  parameters  were  given  by  Yan  and  Alexander  (1990).  These  parameters  were 
chosen  sc  as  to  incorporate  all  known  source  mechanisms  of  the  New  Brunswick 
sequence.  Ideally,  we  would  have  liked  to  compute  synthetics  for  all  possible  ranges 
with  finer  increments,  but  the  training  time  would  have  become  prohibitively  long.  A 
total  of  4608  spectral  ratios  (feature  vectors)  were  calculated  in  the  period  range  of  15  to 
51  sec  (Figure  4). 

Neural  Network  Simulation 

Following  Patnaik  (1990),  several  feed-forward  networks  with  1,  2,  or  3  hidden 
layers  and  varying  numbers  of  processing  units  in  each  of  these  layers  were  simulated 
during  our  experiments.  The  input  layer  in  each  of  these  networks  feeds  a  spectral  ratio 
as  an  input  vector  and  the  8  output  units  correspond  to  the  8  focal  depths.  An  example 
3-layer  (Input-Hidden-Output)  network  simulated  with  a  subset  of  the  set  of  4608 
feature  vectors  is  shown  in  Figure  5.  Note  the  correct  classification  of  all  8  classes 
corresponding  to  the  8  depths  (3  km  to  17  km).  The  final  design  of  the  network,  which 
has  been  used  for  subsequent  focal  depth  estimation,  has  4  layers  consisting  of  12,  60, 
48,  and  8  units.  During  the  training  at  intermediate  steps,  the  ability  of  this  network  to 
correctly  classify  training  feature  vectors  was  repeatedly  tested.  The  results  are  shown 
in  Figure  6.  Note  that  within  about  700  cycles  (where  all  the  patterns  have  been  fed- 
forward  once,  errors  propagated  backward  once,  and  connections  adjusted  to,  consti¬ 
tutes  a  cycle),  the  network  learned  to  classify  correctly  at  nearly  the  77%  level.  As  the 
training  continued,  the  total  mean-square-errors  computed  across  all  the  8  output  units 
decreased,  as  shown  in  Figure  7.  This  curve,  also  known  as  network  learning  curve, 
implied  that  the  gradient  still  existed  and  that  there  had  been  no  encounter  with  local 
minima.  This  meant  that  with  further  training  the  accuracy  or  percent  correct  classifica¬ 
tion  would  increase.  For  these  700  cycles  of  learning,  with  the  back-error-propagation 
algorithm  (1-4),  the  time  spent  were  nearly  23  hours  of  CPU-time  on  a  SUN-4  Sparcsta- 
tion  (12  MIPS).  On  another  simulation  of  a  3-layer  network  vrith  the  mean-field-theory 
learning  algorithm  (5-13),  the  training  process  seemed  fairly  slow.  Since  we  did  not 
encounter  any  local  minima  in  the  back-error-propagation  training,  the  simulation  was 
continued  vrith  the  same  4-layer  network.  When  the  entire  training  set  (4608  features) 
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was  tested  at  the  80%  level,  the  response  of  the  output  units  showed  correct  classifica¬ 
tions  as  shown  in  Figure  8.  The  response  curves  peak  at  the  corresponding  depth 
values. 

The  New  Brunswick  Sequence 

The  January  9,  1982  New  Brunswick  main-shock  and  three  of  its  aftershocks  stu¬ 
died  in  this  paper  are  listed  in  Table  1.  The  locations  of  these  events  and  their  paths  to 
GDSN  stations  are  shown  in  Figure  9.  The  main  shock  source  parameters  have  been 
studied  by  many  investigators  using  body  wave  and  surface  wave  data  {e.g. ,  Choy 
et  al. ,  1983;  Nguyen,  1985;  Wetmiller,  et  at. ,  1984;  Yan  and  Alexander,  1990).  The 
estimated  depth  range  varies  from  7  to  9.5  km  in  these  studies. 

Long-period  surface  wave  time  signals  were  obtained  from  the  GDSN  stations 
ANMO,  ALQ,  SCP,  and  LON.  The  instrument  responses  were  removed  in  the  fre¬ 
quency  domain.  Rayleigh  waves  were  extracted  from  the  vertical  components  and  the 
horizontal  components  were  rotated  to  obtain  Love  waves  (transverse  component). 
Spectral  ratios  between  15  -  51  sec  were  constructed  using  the  phase-matched-filtered 
fundamental-mode  Love  and  Rayleigh  wave  amplitude  spectra.  Upon  testing  with  the 
synthetically  trained  neural  network,  the  estimated  depth  came  out  to  be  7±1  km  for 
the  main  shock.  Figure  10  shows  the  output  response  obtained  from  the  observed  spec¬ 
tral  ratios.  The  peak  corresponds  to  a  target  depth  of  7  km  and  the  broadened  pulse 
shape  provides  the  uncertainty,  which  is  ±1  km  in  this  case. 

Note  that  the  small  hump  (Figure  10)  to  the  left  of  the  peak  at  7  km  is  a  nominal 
response  at  3  km  produced  by  the  spectral  ratio  observed  at  the  station  LON.  We  attri¬ 
bute  this  to  noise,  since  the  majority  of  the  best-matches  peak  at  7  km. 

We  also  conducted  a  least-squares  search  to  find  the  theoretical  ratio  which 
matches  the  main  shock  spectral  ratio  at  station  ANMO.  The  best  fit  is  shown  in  Figure 
11.  The  depth  corresponding  to  this  theoretical  spectral  ratio  is  7  km.  Instead  of  an 
exhaustive  search,  the  search  range  was  limited  to  within  ±10®  of  the  known  main- 
shock  fault  parameters.  The  spectral  ratios  corresponding  to  this  search  range  also  con¬ 
stitutes  a  subset  of  the  original  training  set  of  4608  spectral  ratios. 

Similar  procedures  for  signal  processing  were  adopted  for  the  three  aftershocks 
(Table  1)  of  the  New  Brunswick  event,  and  the  corresponding  spectral  ratios  were  gen¬ 
erated.  Testing  was  done  with  the  same  trained  4-layer  neural  network.  The  output 
responses  obtained  are  also  displayed  in  Figure  10.  The  depths  estimated  for  the  three 
aftershocks  agree  within  ±1  km  of  those  reported  by  Yan  and  Alexander  (1990).  A  com¬ 
parative  listing  of  our  estimated  depths  is  given  in  Table  2. 

Southern  Indian  Ocean  Events 

Since  the  estimated  depths  for  the  New  Brunswick  sequence  were  within  ±1  km  of 
those  reported  by  several  authors,  the  applicability  of  the  present  technique  seems  satis¬ 
factory.  Therefore,  our  proposed  algorithm  (Figure  2)  for  depth  estimation  with 
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multilayered  neural  networks  may  be  utilized  for  events  that  have  not  been  studied.  To 
examine  this  assertion,  we  estimated  the  depths  of  two  oceanic  events  which  occurred 
in  Southern  Indian  Ocean  (Figure  14).  Using  an  earlier  derived  shear  velocity  model 
from  the  inversion  of  Love  and  Rayleigh  wave  group  velocities  and  using  the  same  tech¬ 
niques  for  the  computation  of  Love  to  Rayleigh  wave  spectral  ratios  as  in  the  case  of 
New  Brunswick  sequence,  3240  spectral  ratios,  within  the  period  range  of  18  -  78  sec, 
were  computed  (Figure  12).  These  spectral  ratios  correspond  to  five  different  focal 
depths  (5.0,  7.5,  10.0,  12.5,  and  15  km).  This  chosen  depth  range  is  typical  of  the  focal 
depths  of  events  that  occur  jn  this  region.  The  fault  parameters  chosen  were,  dip  from 
0®  to  80®,  at  10®  intervals,  slip  from  40®  to  150®,  at  10®  intervals  and  strike  from 
175®  to  225®,  at  10®  intervals.  With  the  s)mthetic  spectral  ratios,  several  neural  networks 
with  various  architectures  (combinations  of  number  of  layers  and  number  of  processors 
in  each  layer)  were  simulated.  Irutial  training  with  the  mean-field-theory  learning  algo¬ 
rithm  turned  out  to  be  slow,  as  expected.  Finally,  a  4-layer  feed-forward  back-error- 
propagation  network  (15,  60,  30,  and  5  units)  was  chosen  for  the  training.  The  network 
learning  ctirve  after  5200  cycles  is  shown  in  Figure  13.  The  time  required  for  this  partial 
training  (at  the  90%  level)  was  nearly  37  hours  of  CPU  time  on  a  SUN-4  Sparcstation.  As 
seen  in  Figure  13,  the  gradient  still  exists  and  the  network  could  learn  at  higher  than  the 
90%  level  with  more  training.  The  same  procedure  followed  in  the  previous  simulation 
for  training  and  testing  were  also  followed  in  the  present  experiment. 

The  two  events  chosen  for  depth  estimation  are  shown  in  Figure  14  and  listed  in 
Table  3.  Long  period  SRO  seismograms  were  retrieved  from  the  GDSN  station  NWAO. 
These  were  processed,  and  Love  and  Rayleigh  wave  spectral  ratios  were  generated. 
After  testing  with  the  partially  trained  network  for  the  oceanic  model,  the  output 
responses  were  obtained  as  shown  in  Figure  15.  The  estimated  depth  for  the  08  Jun  82 
event  is  10±1  km  which  agrees  well  with  a  depth  of  10  km  reported  both  by  NEIS  and 
ISC.  The  estimated  depth  for  the  other  event  of  10  Apr  82  is  7.5±1  km  which  is 
between  the  Harvard  Centroid  depth  (10  km)  and  ISC  quoted  depth  5±1.7km,  deter¬ 
mined  from  the  time  difference  of  pP  -  P.  In  general,  the  agreements  are  quite  satisfac¬ 
tory  (Table  4). 

Noise  Tolerance  Test 

Since  the  synthetic  spectral  ratios  used  for  network  training  were  noise  free,  we 
examined  the  robustness  of  the  learning  technique  in  the  presence  of  noisy  feature  vec¬ 
tors  (spectral  ratios).  Using  the  same  oceanic  velocity  model,  300  Love  to  Rayleigh  wave 
spectral  ratios  were  generated  for  three  focal  depths:  5,10,  and  15  km,  corresponding  to 
the  fault  parameters:  dip  10®  to  40®,  at  10®  intervals,  slip  10®  to  50°,  at  10°  intervals,  and 
strike  10®  to  50®,  at  10®  intervals.  The  obvious  reason  for  the  small  number  of  patterns 
is  to  construct  a  smaller  network,  and  thereby,  speed  up  the  training  time.  A  three-layer 
feed-forward  (15,  45,  and  3  units)  network  was  designed  and  trained  with  this  noise-free 
data  set.  Pseudo-random  noise  was  generated  at  10%,  20%,  j07o  and  40%  levels  and 
added  to  the  spectral  ratios.  Each  of  the  four  sets  of  these  noise-contaminated  spectral 
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ratios  were  then  used  to  train  and  test  the  same  network.  The  network  learning  curves 
for  all  five  cases  are  shown  in  Figure  16.  The  solid  line  is  the  noise-free  data,  and  open 
circles,  triangles,  pluses,  and  crosses  correspond  to  10%,  20%  30%  and  40%  noise  lev¬ 
els,  resf)€ctively.  As  is  evident  from  Figure  16,  the  total  mean-square-error  remains 
higher  for  the  training  data  with  higher  noise  level.  However,  the  important  observa¬ 
tion  is  that  even  with  a  40%  noise  level,  the  network  still  learns  to  correctly  classify  over 
70%  of  the  training  set.  The  gradient  in  each  of  the  curves  indicates  that  with  further 
training,  better  performance  could  be  obtained.  Figure  17  shows  the  histograms  of  the 
percentage  of  correct  classifications  for  each  target  depth.  The  histograms  show  results 
obtained  for  testing  the  training  sets  with  different  noise  levels,  after  2000  cycles  each. 


Discussion  and  Conclusions 

This  study  has  shown  that  neural  networks  can  successfully  be  used  with  surface 
wave  spectral  ratios  to  obtain  focal  depths  of  shallow  events.  During  our  initial  simula¬ 
tion  and  subsequent  testing  with  observed  data  it  was  seen  that  success  of  the  resultant 
networks  required  accurate  source  and  crustal  models.  This  was  also  evident  in  earlier 
studies  directly  using  spectral  ratios  {e.g. ,  Harkrider,  1970;  Masse  et  al. ,  1973;  Van  and 
Alexander,  1990).  The  neural  network  technique  developed  here  uses  the  location  of 
the  spectral  holes  {e.g. ,  Harkrider,  1970;  Tsai  and  Aki,  1970;  Masse  et  al. ,  1973)  as  well 
as  the  shape  of  the  normalized  amplitude  spectral  ratios  to  obtain  a  map  between 
depths  and  corresponding  ratios.  Azimuth  dependence  is  included  through  the  choice 
of  many  strike  possibilities  of  the  source.  Although  all  static  source  parameters  (depth, 
dip,  slip,  strike)  can,  in  principle,  be  obtained  with  this  technique,  because  of  the  CPU¬ 
time  requirements,  we  restricted  our  problem  only  to  depth  estimation.  Determination 
of  all  parameters  will  require  a  faster  computer  such  as  a  parallel  processor  to  speed  up 
the  training  process.  Once  trained,  however,  the  network  can  be  implemented  on  a 
slower  machine  for  application  purposes. 

One  problem  with  this  technique  is  that,  at  present,  there  is  no  network  architec¬ 
ture  available  which  can  account  for  missing  data  points;  therefore,  interpolated  values 
must  be  supplied  to  fill  in  those  points.  It  is  therefore  desirable,  when  making  determi¬ 
nations  of  focal  depth,  to  use  observed  spectral  frequencies  which  are  identical  to  those 
used  for  training,  so  that  assumed  values  do  not  need  to  be  provided. 

In  the  noise  tolerance  test  we  have  added  random  noise  to  the  surface  wave  spec¬ 
tral  data.  For  deterministic  signals  (well-dispersed  surface  waves)  such  noise  may  not  be 
realistic.  Therefore  an  alternative  approach  might  be  to  obtain  actual  noise  from  long 
period  digital  seismograms,  over  the  time  interval  preceding  the  signal,  and  add  that  to 
the  synthetic  data  for  training. 

The  generalization  capability  of  a  trained  neural  network  for  the  problem  of  focal 
depth  estimation  must  be  addressed,  where  generalization  means  reasonable  response 
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to  unknown  inputs.  Such  inputs  may  fall  outside  the  mapping  boundary  of  the  network. 
As  an  example,  one  of  the  networks  has  been  trained  for  the  depth  range  of  3  km  to  17 
km.  Therefore,  the  network  is  capable  of  approximating  those  spectral  features  that  are 
continuously  mapped  in  this  depth  range.  If,  however,  a  spectrcd  ratio  corresponding  to 
a  depth  of  30  km  is  presented  to  the  network,  it  would  produce  some  response  which 
would  be  meaningless.  This  difficulty  arises  from  the  fundamental  problem  inherent  in 
mappings  through  bounded  inputs/outputs.  Since  the  feature  space  (possible  range  of 
focal  depths  in  a  region)  is  unbounded  and  the  output  classes  are  bounded  (3  km  to  17 
km),  these  difficulties  may  arise.  However,  there  can  always  be  judicious  applications  of 
this  useful  technique.  From  a  seismological  perspective,  one  may  first  choose  to  use 
deeper  event  depth  identifiers  such  as  pP-P  times.  If  that  depth  estimate  falls  outside 
the  depth  range  of  a  trained  network  it  should  not  be  used.  The  present  technique  is 
suitable  for  shallower  events  only.  On  the  other  hand,  one  may  choose  to  synthesize 
features  for  a  larger  depth  range  to  train  a  network.  This  will  be  limited  again  by  the 
availability  of  computing  power. 

Currently  there  is  no  well-established  way  to  predict  the  internal  dynamics  of  a 
neural  network  so  as  to  understand  how  a  network  actually  learns  a  given  mapping. 
Approximate  techniques  such  as  the  analysis  of  the  connection  matrix  of  a  trained  net¬ 
work  is  often  done  for  simpler  problems  and  smaller  networks  (Elman  and  Zipser,  1988; 
Gorman  and  Sejrwwski,  1988).  For  the  present  application  of  a  4-layer  network,  such  an 
analysis  would  be  extremely  difficult,  and  perhaps  meaningless.  From  a  physical  per¬ 
spective,  the  network  dynamics  is  seen  as  chaotic  and  attempts  have  been  made  (e.g. , 
Choi  and  Huberman,  1983;  Geman  and  Geman,  1984;  Clark,  1988)  to  study  the  underly- 
ing  physics.  The  theory  of  neural  computing  is  still  evolving,  and  many  issues  still 
remain  to  be  answered. 

In  this  study  involving  neurcil  computing,  we  have  presented  a  data-intensive, 
region-specific  novel  technique  of  nonlinear  transfer  function  approximation  as  applied 
to  the  estimation  of  earthquake  focal  depths  from  surface  wave  amplitude  spectral 
ratios.  The  technique  is  to  train  a  neural  network  with  the  depth-dependent  mapping  of 
the  spectral  ratios  such  that  when  an  observed  spectral  ratio  is  processed  through  it,  a 
corresponding  depth  is  estimated. 

Application  of  the  technique  to  the  1982  New  Brunswick  main-shock  and  three  of 
its  aftershocks  produced  depth  estimates  within  ±1  km  of  the  well  calibrated  depths 
reported  by  other  investigators. 

Application  to  two  events  in  the  southern  Indian  Ocean  provided  depths  for 
which,  one  of  the  events  (04/10/82)  agreed  to  within  ±2.5  km  of  the  reported  centroid 
depth  and  ISC  reported  depth  from  the  (pP  —  P)  times.  For  the  other  event  (06/08/82), 
the  estimated  depth  was  within  ±1  km  of  the  depths  reported  for  that  region  by  the 
NEIS  and  ISC. 
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Figure  Captions 


Fig  1.  The  schematic  on  the  left  represents  a  three-layer  neural  network  which  can  be 
trained  using  a  learning  algorithm  (Error-back-propagation,  Boltzmaim 
Machine,  or  Mean-field-theory).  TT\e  upper-right  figure  shows  an  analog 
equivalent  of  the  functions  within  a  single  neuron  (unit).  The  non-linear  ele¬ 
ment  uses  a  sigmoidal  function  (bottom-right).  Note  that  at  0  activation  and  high 
gain,  the  neuron  behaves  as  a  comparator  or  a  flip/fiop.  However,  at  low  gain 
and  other  activation  values,  the  non-linearity  controls  the  overall  dynamics  of 
the  network. 

Fig  2.  A  schematic  representation  of  focal  depth  estimation  with  a  neural  network.  The 
neural  network  is  trained  with  synthetic  data  and  tested  with  observed  and/or 
synthetic  data. 

Fig  3.  The  shear-velocity  model  of  Taylor  and  Toksoz  (1982)  for  the  Eastern  United 
States.  This  model  is  used  for  the  calculation  of  synthetic  surface  wave  spectral 
ratios. 

Fig  4.  Superposition  of  4608  normalized  spectral  ratios  calculated  with  the  Eastern 
United  States  velocity  model  for  different  combinations  of  dip,  slip,  strike  and 
depth. 

Fig  5.  Schematic  of  an  example  3-layer  feed-forward  back-error-propagation  neural  net¬ 
work  that  has  been  trained  with  a  shbset  of  Love  to  Rayleigh  wave  spectral 
ratios  computed  for  the  eastern  United  States  velocity  model.  The  connection 
matrices  are  in  the  middle.  The  Input,  Hidden  and  Output  responses  are 
shown  at  the  right. 

Fig  6.  Improvements  in  network  learning  by  increasing  the  number  of  training  steps. 
The  histograms  show  the  number  of  patterns  classified  at  a  given  target  depth. 
The  darker  bars  are  the  final  network  estimates.  The  lighter  shaded  bars  are  the 
results  obtained  from  intermediate  tests  using  50  and  450  steps.  Note  the  wider 
distribution  of  depths  at  early  learning  stages  of  the  network.  When  the 
number  of  training  steps  increases,  the  depth  distribution  for  corresponding 
target  depths  become  sharper.  The  expected  distribution  for  a  classification  that 
is  100%  correct  is  a  spike  located  at  that  target  depth. 

Fig  7.  The  network  learning  curve  (total  mean-square-error  vs.  number  of  training 
cycles)  of  the  neural  network  trained  with  the  spectral  ratios  computed  from  the 
eastern  United  States  velocity  model. 
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Fig  8.  Response  of  the  partially  trained  network  upon  testing  with  the  training  set.  Each 
response  curve  peaks  at  its  corresponding  depth.  The  sharper  the  peaks,  the 
better  the  resolution. 

Fig  9.  Location  of  the  New  Brunswick  earthquake  and  its  aftershocks.  The  DWWSSN 
stations  are  LON,  ALQ,  and  SCP.  ANMO  is  the  SRO  station. 

Fig  10.  Response  of  the  partially  trained  network  (over  80%  learned)  to  the  spectral 
ratios  of  the  New  Brunswick  earthquake  and  three  of  its  aftershocks.  The 
peaks  correspond  to  estimated  depths. 

Fig  11.  A  comparison  of  normalized  observed  (solid  circles)  and  theoretical  (solid  line) 
fundamental  mode  amplitude  spectral  ratio  of  Love  to  Rayleigh  wave  for  the 
Eastern  United  States  model. 

Fig  12.  Superposition  of  3240  normalized  spectral  ratios  calculated  for  the  southeastern 
Indian  Ocean  velocity  model  for  different  combinations  of  dip,  slip,  strike  and 
depth. 

Fig  13.  Network  learning  cxurve  corresponding  to  the  southern  Indian  Ocean  velocity 
model.  The  high-frequency  oscillations  on  this  curve  are,  possibly,  the  result  of 
underlying  dynamics  which  is  often  chaotic. 

Fig  14.  Epicenters  of  the  events  located  in  the  southern  Indian  Ocean  and  their  paths  to 
the  GDSN  station  NWAO. 

Fig  15.  Response  of  the  partially  trained  network  (over  90%  learned)  to  the  spectral 
ratios  of  two  events  in  the  southern  Indian  Ocean.  The  peaks  correspond  to  the 
estimated  depths. 

Fig  16.  Network  learning  curve  with  respect  to  training  at  various  levels  of  noise  added 
to  the  spectral  ratios. 

Fig  17.  The  nobe  tolerance  test  conducted  with  a  smaller  network  for  5,10,  and  15  km 
target  depths.  The  histograms  show  the  depth  distribution  over  each  target 
depth  with  10%,  20%,  30%,  and  40%  noise  levels.  The  darkest  bars  represent 
the  10%  noise  level  and  the  lightest  bars  represent  the  40%  noise  level.  When 
noise  level  is  increased,  the  percentage  of  correct  classifications  gradually 
decreases.  At  all  the  noise  levels,  the  percentage  remains  almost  the  same 
which  indicates  that  the  network  is  relatively  noise  insensitive. 
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Figure  17 


Chapter  6 


Discrimination  Studies  Using  Regional 
Phases  in  the  Central  United  States 


by 


M.  Amjad  and  BJ.  Mitchell 


ABSTRACT 

The  usefulness  of  amplitude  ratios  discriminants  in  the  central  United 
States  is  systematically  studied  using  signals  from  small  earthquakes  and 
explosions  recorded  at  two  broad-band  stations  operated  by  Saint  Louis 
University.  Signals  are  partitioned  into  passbands  of  0-2  Hz,  2-4  Hz,  4-6  Hz, 
and  8-10  Hz.  Ratios  of  different  phases  to  one  another  and  of  one  passband 
to  another  for  the  same  phase  were  tested  for  75  events  presumed  to  be 
explosions  and  40  events  presumed  to  be  earthquakes  recorded  at  station 
FVM  over  the  period  May  983  -  June  1989.  Ratios  were  determined  for  the 
whole  signal  as  well  as  for  the  phases  Pg,  Sg,  Lg,  and  Rg.  Some  of  the  ratios 
appear  to  provide  useful  discriminants,  with  the  ability  to  discriminate 
increasing  with  increasing  frequency. 

Before  correcting  for  instrument  response,  the  amplitudes  for  earth¬ 
quakes  for  the  whole  signal,  Lg,  and  Rg  tend  to  increase  over  the  range  1-9 
Hz  and  tend  to  decrease  for  explosions  over  the  same  range.  None  of  the  75 
explosions  were  found  to  have  increasing  amplitudes  from  0  to  9  Hz  on  all 
three  components  for  all  three  phases.  Only  4  out  of  40  earthquakes  were 
found  to  have  a  decreasing  trend  on  all  three  components  for  amplitudes  of 
the  same  three  phases.  This  variation  of  amplitudes  is  not  affected  by  dis¬ 
tance  at  distances  between  20  and  340  km,  or  by  source  strength. 

6  known  explosions  and  5  known  earthquakes  were  recorded  at  the  IRIS 
station  CCM  over  the  period  September  1989  to  November  1990.  Using  the 
criteria  for  discrimination  suggested  by  the  FVM  data,  all  of  these  events 
would  be  correctly  identified. 

After  correction  for  instrument  responses,  ratios  at  both  stations  FVM 
and  CCM  exhibit  decreasing  values  of  ground  motion  with  increasing  fre¬ 
quency,  but  the  slope  for  explosions  is  significantly  greater  than  that  for 
explosions. 


INTRODUCTION 

For  the  past  few  decades  the  focus  of  seismic  research  for  verification  of 
nuclear  test  ban  treaties  has  been  to  develop  methods  for  estimating  the  yield 
of  nuclear  explosions  from  seismic  records  at  large  distances.  Seismic  moni¬ 
toring  networks  were  installed  to  enhance  the  capabilities  for  better  recording 
at  large  distances.  Research  related  to  nuclear  test-ban  treaties  led  to  a 
number  of  discriminants  for  out-country  verificiation  (Pomeroy  et  al.,  1982), 
and  the  methods  for  estimating  yields  of  large  nuclear  explosions  (Blandford, 
1982). 

Recently,  there  have  been  new  developments  in  test  ban  treaty  negotia¬ 
tions  between  the  U.S.S.R.  and  the  U.S.A.  These  new  developments  in 
negotiations  demand  on-site  verification,  and  deployment  of  in-country 
seismic  monitoring  capabilities  for  verification  from  in-country  seismi, 
records.  In  response  to  these  developments  in  negotiations,  and  in  anticipa¬ 
tion  of  a  range  of  new  test  ban  treaties,  nuclear  test  ban  research  has  taken 
new  directions.  These  include  attempts  to  find  discriminants  between  small 
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nuclear  explosions  and  earthquakes,  and  to  develop  methods  for  estimating 
yield  and  source  parameters  of  low  yield  events.  Further  reserach  is  likely  to 
utilize  suitable  instruments  at  short  distances,  and  will  focus  on  finding 
appropriate  discriminants  between  small  explosions  and  small  earthquakes, 
and  on  finding  methods  for  estimating  yield  of  small  nuclear  explosions. 
This  study  attempts  to  find  suitable  discrimination  criteria  to  separate  chemi¬ 
cal  explosions  and  earthquakes  in  the  central  United  States.  This  region  has 
been  found  to  be  similar  to  large  areas  of  the  Soviet  Union  for  wave  propaga¬ 
tion  studies  (Sweeney  et  al.,  1986). 

The  New  Madrid  seismic  zone  is  a  seismically  active  region  of  the  central 
United  States.  Most  of  the  seismic  activity  occurs  in  a  small  region;  however, 
many  earthquakes  occur  in  a  broader  region  surrounding  the  New  Madrid 
fault  system.  There  is  also  much  mining,  and  quarrying  activity  taking  place 
in  the  same  region.  Most  of  the  mining  explosions  and  quarry  blasts  occur  as 
strip  mining  for  coal  and  quarrying  for  coal  and  quanying  for  limestones,  but 
other  forms  of  surface  mining  also  occur.  In  addition,  underground  mining 
activities  for  iron,  lead,  and  coal  are  conducted. 

Figure  1  shows  the  locations  of  earthquakes  and  explosions  utilized  in 
this  study.  They  occur  in  the  states  of  Missouri,  Illinois,  Indiana  and  Ken¬ 
tucky  and  were  recorded  in  digital  form  at  two  stations,  FVM  and  CCM, 
operated  by  Saint  Louis  University.  The  size  of  the  explosions  for  quarry 
blasts  range  between  1  ton  and  10  tons,  and  the  size  of  the  strip  mine  explo¬ 
sions  range  between  1  ton  and  50  tons,  with  the  exception  of  one  or  two 
explosions  of  about  100  Ktons.  10  to  30  small  explosions  occur  every  day. 
Data  is  available  covering  the  distance  range  5-340  km,  both  for  earthquakes 
and  chemical  explosions.  Comparative  studies  of  signal  characteristics  can 
therefore  be  conducted  using  these  sources.  Station  FVM  has  operated  in  a 
digital  mode  since  1980.  CCM  is  an  IRIS  station  (Incorpprated  Research  Insti¬ 
tutions  for  Seismology)  installed  by  Saint  Louis  University  in  1990. 

The  character  of  recorded  ground  motion  depends  on  source  parameters, 
instrument  response,  and  the  physical  properties  of  the  medium  through 
which  seismic  waves  propagate.  The  physical  properties  of  the  earth  affect 
the  seismic  signal  in  various  ways,  e.g.  by  anelastic  attenuation,  scattering, 
and  velocity  structure  of  the  crust.  These  physical  properties  are  determined 
ultimately  by  the  geology  and  tectonic  setting  of  the  region.  Hence  seismic 
records  obtained  in  different  regions  differ  in  signal  characteristics.  It  has 
been  realized  that  discriminants  which  are  useful  in  one  region  may  not  work 
in  others. 

The  results  of  previous  studies  indicate  that  regional  discriminants  can¬ 
not  be  used  with  the  same  effectiveness  everywhere.  There  seems  to  be  a 
strong  dependence  on  near-source  geology  and  propagation  path,  which  vary 
from  region  to  region.  Because  of  the  strong  regional  bias  in  discrimination 
criteria,  there  is  a  need  for  a  case-based  approach  to  identify  an  unknown 
event,  whethpr  it  is  an  explosion  or  earthquake.  The  most  useful  approach 
would  be  to  compare  the  characteristics  of  an  unknown  event  with  the 
characteristics  of  known  events  in  the  same  region.  There  is,  therefore,  a 
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need  to  find  viable  discriminants  using  known  events  from  different  geologi¬ 
cal  regions.  Baumgardt  and  Young  (1990)  have  discussed  the  importance  of  a 
case  based  approach  for  discrimination. 

Data  Collection  and  Processing 

Data  were  collected  for  events  at  epicentral  distances  between  5  and  340 
km  from  FVM  and  between  134  and  417  km  fi-om  CCM.  Seismograms  were 
obtained  for  75  presumed  explosions  and  40  presumed  earthquakes  recorded 
by  the  station  FVM  from  1983  to  1988  and  by  5  known  earthquakes  and  6 
known  explosions  recorded  at  CCM  in  1989  and  1990.  The  locations  of  the 
stations  and  events  are  shown  in  Figure  1.  Since  there  is  no  explosion  data 
from  the  seismically  active  New  Madrid  seismic  zone,  no  earthquakes  were 
selected  from  that  region.  The  magnitudes  of  the  explosions  range  from  1.9 
to  3.4  and  the  magnitudes  of  the  earthquakes  from  1.1  to  3.9.  Example 
vertical-component  records  at  station  FVM  for  some  of  the  events  appear  in 
Figure  2. 

We  rotated  the  horizontal  component  records  into  radial  and  transverse 
components.  The  onsets  of  the  Pg,  Sg,  Pn  and  Sn  phases  were  then  picked 
manually  and  a  one-second  window  was  used  from  the  start  of  the  phase 
selected  for  analysis.  This  window  was  applied  to  avoid  contamination  from 
closely  following  phases.  For  the  Lg  phase,  a  fixed  velocity  window  from 
3.60  km/sec  to  3.00  km/sec  was  used.  For  each  event,  radial,  transverse  and 
vertical  component  signals  were  filtered  to  obtain  records  in  the  frequency 
ranges  0  to  2  Hz,  2  to  4  Hz,  4  to  6  Hz,  6  to  8  Hz,  and  8  to  10  Hz,  using  a 
zero-phase,  two-pole,  time-domain  Butterworth  filter. 

Tests  of  Possible  Discriminants 

We  looked  at  three  possible  discriminants:  (1)  amplitude  ratios  using 
various  phases,  (2)  ratios  between  lower  and  higher  frequencies  utilizing  the 
same  phase,  and  (3)  the  variation  of  amplitudes  with  frequency  for  individual 
seismograms.  A  summary  of  the  significant  results  is  given  in  the  following 
paragraphs.  A  complete  presentation  of  those  results  appears  in  Amjad 
(1991). 

Amplitude  ratios  of  various  phases  have  been  studied  by  several  authors 
(e.g.  Nuttli,  1981;  Pomeroy  et  al.,  1982;  Bennett  and  Murphy,  1986;  Taylor  et 
al.,  1989;  Bennett  et  al.,  1989)  with  mixed  results.  In  the  present  study,  the 
ratios  between  different  phases  over  the  same  frequency  passband  for  the 
same  component  are  calculated  for  Pg/Sg,  Pg/Lg,  Pg/Pn,  Pg/Rg,  Sg/Lg,  Sg/Rg, 
Sg/Pn,  Lg/Rg,  Lg/Pn,  and  Rg/Pn  phases.  There  are  10  possible  ratios  of  dif¬ 
ferent  phases  e.g.  Pg/Sg  for  each  frequency  range.  All  these  possible  ratios 
for  radial,  transverse,  and  vertical  components  were  determined  as  functions 
of  distance  and  magnitude. 

Although  there  is  some  suggestion  of  separation  between  earthquakes 
and  explosions  for  some  ratios  in  certain  passbands,  the  separation  is  never 
clear-cut.  The  number  of  separated  events  vary  from  component  to  com¬ 
ponent  for  different  phase  ratios  and  for  different  frequency  ranges.  For 
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instance,  the  ratio  Sg/Rg  indicates  some  separation  at  higher  frequencies, 
whereas  the  ratio  Pg/Rg  inhibits  better  separation  at  low  frequencies.  In  most 
cases,  if  there  is  separation,  it  tends  to  increase  for  higher  frequency 
passbands  of  6-8  Hz  and  8-10  Hz.  Our  results  generally  indicate,  however, 
that  the  ratio  of  different  phases  with  respect  to  each  other  is  not  an  effective 
discrimination  criterion  in  the  central  United  States. 

Ratios  between  lower  and  higher  frequencies  for  regional  phases  were 
studied  by  Murphy  and  Bennett  (1987)  and  Taylor  et  al.  (1988)  and  were 
found  to  be  useful  in  many  cases.  Bennett  et  al.  (1989),  however,  reported 
less  effectiveness  in  some  regions  than  in  others.  In  this  study,  maximum 
amplitudes,  rms  values,  maximum  spectral  amplitudes,  mean  spectral  ampli¬ 
tudes,  and  rms  spectral  Vcdues  of  Pg,  Sg,  Pn  phases,  and  values  of  maximum 
amplitudes  of  whole  signal,  of  Lg  and  Rg  phases  were  measured.  These 
measurements  are  made  for  the  whole  signal  (0-10  Hz),  and  for  0-2  hz,  2-4 
Hz,  4-6  Hz,  6-8  Hz,  and  8-10  Hz  frequency  passbands  for  radial,  transverse 
and  vertical  components  for  all  events.  The  0-10  Hz/0-2  Hz,  0-10  Hz/2-4  Hz, 
0-10  Hz/4-6  Hz,  0-10  Hz/6-8  Hz,  0-10  Hz/8-10  Hz,  0-2  Hz/2-4  Hz,  0-2  Hz/4-6 
Hz,  0-2  Hz,  6-8  Hz,  0-2  Hz/8-10  Hz,  2-4  Hz/4-6  Hz,  2-4  hz/6-8  Hz,  2-4  Hz/8-10 
Hz,  4-6  Hz/6-8  Hz,  4-6  Hz/8  10  Hz,  and  6-8  Hz/8-10  Hz  ratios  for  Pg,  Sg,  Lg, 
Rg,  and  Pn  were  calculated  for  radial,  transverse  and  vertical  components. 
These  ratios  indicated  some  separation  between  earthquakes  and  explosions 
could  be  obtained,  with  best  results  occurring  for  the  Lg  and  Rg  phases. 

Figure  3  illustrates  results  for  the  Rg  phase  recorded  at  FVM.  The  low 
frequency  ratio  places  most  explosions  below  most  earthquakes  whereas  the 
high-frequency  ratio  shows  the  opposite  relation.  There  is  some  overlap 
which  cannot  be  explained,  but  it  is  possible  that  some  events  have  been 
misidentified.  A  difficulty  with  Rg  is  that,  while  it  can  be  a  large  phase  for 
explosions,  it  can  also  be  readily  absorbed  on  transmission  through  low-Q 
material.  This  seems  to  have  occurred  for  several  of  the  explosion  events 
shown  in  Figure  2.  Similar  separations  in  ratios  in  Figure  3  were,  however, 
also  observed  for  Lg  as  recorded  on  all  three  components. 

The  most  promising  discriminant  which  we  found  relies  on  the  variation 
of  amplitudes  with  frequency.  After  correcting  for  instrument  response,  the 
amplitude  decreases  with  increasing  frequency  is  always  greater  than  that  for 
explosions.  Figure  4  shows  results  for  Lg  waves  using  vertical-component 
recordings.  Similar  results  were  found  for  radial  and  transverse  components 
of  Lg  and  for  all  components  of  Rg. 

Conclusions 

Using  ratios  of  different  passbands  for  the  same  phase  is  more  promising 
than  using  ratios  of  different  phases  for  the  same  passband  in  discriminating 
between  earthquakes  and  explosions  in  the  central  United  States.  The  varia¬ 
tion  of  amplitudes  with  frequency  between  the  0-2  Hz  and  6-8  Hz  passbands, 
however,  is  the  most  promising  discriminant  for  found  in  our  study.  It 
appears  to  apply  to  all  phases,  with  best  results  likely  to  be  obtained  for  the 
Lg  and  Rg  phases. 
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Figure  Captions 


Figure  1.  Map  of  area  of  study  showing  locations  of  seismic  stations  (FVM 
and  CCM),  explosions,  and  earthquakes. 

Figure  2.  Example  vertical-component  seismograms  generated  by  explosions 
(left)  and  earthquakes  (right)  as  recorded  at  station  FVM. 

Figure  3.  Ratios  of  Rg  amplitudes  at  various  frequencies  versus  distance  as 
determined  from  vertical-component  records. 

Figure  4.  Normalized  Lg  amplitude  variations  for  the  vertical  component 
after  correcting  for  instrument  response.  (Top)  Amplitude  variations  at 
CCM  for  explosions  (left)  and  earthquakes  (right).  (Bottom)  Amplitude 
variations  at  FVM  for  explosions  (left)  and  earthquakes  (right). 
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